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Part I 



Single polaron 

I. INTRODUCTION. THE "STANDARD" THEORIES 
A. The polaron concept 

A charge placed in a polarizable medium is screened. Dielectric theory describes the 
phenomenon by the induction of a polarization around the charge carrier. The induced 
polarization will follow the charge carrier when it is moving through the medium. The 
carrier together with the induced polarization is considered as one entity (see Fig.[T]). It was 
called a polaron by L. D. Landau The physical properties of a polaron differ from those 
of a band-carrier. A polaron is characterized by its binding (or self-) energy Eq, an effective 
mass m* and by its characteristic response to external electric and magnetic fields (e. g. dc 
mobility and optical absorption coefficient). 
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FIG. 1: Artist view of a polaron. A conduction electron in an ionic crystal or a polar semiconductor 
repels the negative ions and attracts the positive ions. A self-induced potential arises, which acts 
back on the electron and modifies its physical properties. (From j^.) 

If the spatial extension of a polaron is large compared to the lattice parameter of the solid, 
the latter can be treated as a polarizable continuum. This is the case of a large (Frdhlich) 
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polaron. When the self-induced polarization caused by an electron or hole becomes of the 
order of the lattice parameter, a smll (Holstein) polaron can arise. As distinct from large 
polarons, small polarons are governed by short-range interactions. 



B. Intuitive concepts 

a. The polaron radius. Large polarons vs small polarons Consider the LO phonon field 
with frequency cjlo interacting with an electron. Denote by Au the quadratic mean square 
deviation of the electron velocity. In the electron-phonon interaction is weak, the electron 
can travel a distance 

Ax ^ (1.1) 

during a time Wlq, characteristic for the lattice period, because it is the distance within 
which the electron can be localized using the phonon field as measuring device. From the 
uncertainty relations it follows 



ApAx = [Au) n 



Au '^"^^ 



m ' 



Ax~./-^. (1.2) 

At weak coupling Ax is a measure of the polaron radius Vp. To be consistent, the polaron 
radius Vp must be considerably larger than the lattice parameter a. (this is a criterion of a 
"large polaron"). Experimental evaluation of the polaron radius leads to the follwing typical 
values: rp ^ loA for alkali halides, rp ^ 20^4 for silver halides, rp lOOA for 11- VI, 11-V 
semiconductors. The continuum approximation is not satisfied for transition metal oxides 
(NiO, CaO, MnO), in other oxides (U02,Nb02...). For those solids the "small polaron" 
concept is used. In some substances (e.g. perovskites) some intermediate region between 
large and small polarons is realized. 

b. The coupling constant Consider the case of strong electron-phonon interaction in 
a polar crystal. The electron of mass m is then localized and can - to a first approximation 



6 



- be considered as a static charge distribution within a sphere with radius h. The medium 
is characterized by an average dielectric constant £, which will be defined below. 

The potential energy of a sphere of radius /luniformly charged with the charge e in a 
vacuum is (see Eq. (8.6) of Ref. {4]) 

3p2 

U.ac = (1.3) 

The potential energy of a uniformly charged sphere in a medium with the high-frequency 
dielectric constant £00 is 

This is the potential energy of the self-interaction of the charge e uniformly spread over 
the sphere of radius li in a medium with the dielectric constant Soo- In a medium with 
an inertial polarization field (due to LO phonons), the potential energy of the uniformly 
charged sphere is 

where Eq is the static dielectric constant. The polaron effect is then related to the change of 
the potential energy of the interaction of the charged sphere due to the inertial polarization 
field. This change is the potential energy U2 of the uniformly charged sphere in the presence 
of the inertial polarization field minus the potential energy of the self-interaction f/iof the 
charge e uniformly spread over the sphere in a medium without the inertial polarization: 

^ 5 /i (^^o £^00) 5 ell 

with 

1 _ 1 1 

The electron distribution in a sphere may be non-uniform, what may infiuence the nu- 
merical coefficient in Eqs. (I1.4p to ( II. 6p . In this connection one can use the estimate 

Upol ~ (1.7) 

The restriction of the electron in space requires its de Broglie wave length to be of the 
order /i,so that its kinetic energy is of the order 47r^/),^/2m/p Minimizing the total energy 
with respect to /ileads to 
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wherefrom the binding energy is 



4 

em 



For weak coupling, one can neglect the kinetic energy of the electron. Taking the polaron 
radius according to fll.2p . = \j2h/mui^o, the binding energy is 



We note that 



f/2 = -- = --./^. (1.9) 



Ui e^m 1 / f/r 



^ ' '1.10) 



Following the conventions of the field theory, the self energy at weak coupling is written as 
Therefore the so-called Frohlich polaron coupling constant is 




a 



. . ■ (1-11) 

he y 2hujLo £ 

For the average dielectric constant one shows that 

1 _ 1 1 
_ ^ 

where £00 and Eq are, respectively, the electronic and the static dielectric constant of the polar 
crystal. The difference I/Eoo — I/eq arises because the ionic vibrations occur in the infrared 
spectrum and the electrons in the shells can follow the conduction electron adiabatically. 
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c. Polaron mobility Here we give a simple derivation leading to the gross features of 
the mobility behaviuor, especially its temperature dependence. The key idea is that the 
mobility will change because the number of phonons in the lattice, with which the polaron 

interacts, is changing with temperature. 
The phonon density is given by 

1 

n - 



e kT — 1 



The mobility for large polaron is proportional to the inverse of the number of phonons: 

1 ''"LP 

II ^ — — e — 1 
n 

and for low temperatures kT <C ^lo 

ji^ e kT _ (1-12) 

The mobility of continuum poarons decreases with increasing temperature following an ex- 
ponential law. The slope of the straight hne in In/x vs 1/T is characterized by the LO 
phonon frequesncy. Systematic study performed, in particular, by Frohlich and Kadanoff, 
gives 

= e kT . (1-13) 

The small polaron will jump from ion to ion under the influence of optical phonons. The 
lerger the numver of phonons, the lerger the mobility. The behaviuor of the small polaron 
is the opposite of that of the large polaron. One expects: 



e kT — 1 

For low temperatures kT <^ ^lo one has: 

fjLKi e kT : (1-14) 

the mobility of small polaron is thermally activated. Systematic analysis within the small- 
polaron theory shows that 
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^ e : (1-15) 

with 7 ~ 5. 

C. The Frohhch Hamiltonian 

Frohhch proposed a model Hamiltonian for the "large" polaron through which its dynam- 
ics is treated quantum mechanically ("Frohlich Hamiltonian"). The polarization, carried by 
the longitudinal optical (LO) phonons, is represented by a set of quantum oscillators with 
frequency Wlo^ the long- wavelength LO-phonon frequency, and the interaction between the 



charge and the polarization field is linear in the field js]: 

^ = ^ + E ^^Loa^ak + Y.iV.a^e''^-^ + V^^e"*^-^), (1.16) 
^ k k 

where r is the position coordinate operator of the electron with band mass rrib, p is its 
canonically conjugate momentum operator; and Ok are the creation (and annihilation) 
operators for longitudinal optical phonons of wave vector k and energy hui^o- The Vk are 
Fourier components of the electron-phonon interaction 



The strength of the electron-phonon interaction is expressed by a dimensionless coupling 
constant a, which is defined as: 

^^-i.^ (-'-)■ (1.13) 



In this definition, Eoo and Eq are, respectively, the electronic and the static dielectric constant 
of the polar crystal. 

In Tabled] the Frohlich coupling constant is given for a few solids^. 

In deriving the form of Vk, expressions f ll.l7p and (11.181) . it was assumed that (i) the spatial 
extension of the polaron is large compared to the lattice parameter of the solid ( "continuum" 



^ In some cases, due to lack of reliable experimental data to determine the electron band mass, the values 
of a are not well established. 



10 



approximation), (ii) spin and relativistic effects can be neglected, (iii) the band-electron has 
parabolic dispersion, (iv) in line with the first approximation it is also assumed that the 
LO-phonons of interest for the interaction, are the long-wavelength phonons with constant 
frequency wlo- 

The model, represented by the Hamiltonian fll.l6p (which up to now could not been 



solved exactly) has been the subject of extensive investigations, see, e. g., Refs. j2l|-|28|. 
In what follows the key approaches of the Frohlich-polaron theory are briefly reviewed with 
indication of their relevance for the polaron problems in nanostructures. 



D. Infinite mass model ["shift" operators] 

Here some insight will be given in the type of transformation that might be useful to 
study the Frohlich Hamiltonian fll.l6p . For this purpose the Hamiltonian will be treated for 
a particle with infinite mass m?, — )■ oo, (which is at r = 0) : 

H^ = Y^ ^Lo4«k + J2{Vka^ + V^ai), (1.19) 

k k 

TABLE I: Electron-phonon coupling constants (After Ref. Q]) 



Material 


a 


Ref. 


Material 


a Ref. 


InSb 


0.023 [5] 


AgCl 


1.84 [yj 


InAs 


0.052 [5] 


KI 


2.5 [5] 


GaAs 


0.068 [5] 


TlBr 


2.55 [5] 


GaP 


0.20 


[5] 


KBr 


3.05 [5j 


CdTe 


0.29 


m 


Bii2Si02o 


3.18 [19] 


ZnSe 


0.43 




CdF2 


3.2 [5J 


CdS 


0.53 




KCl 


3.44 [5J 


a-Al203 


1.25 


m 


Csl 


3.67 [5J 


AgBr 


1.53 


m 


SrTiOg 


3.77 [20] 


a-Si02 


1.59 


[12] 


RbCl 


3.81 [5] 
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which can be transformed into the following expression with "shifted" phonon operators: 

|2 



v: 



-E 



(1.20) 



To determine the eigenstates of this Hamiltonian, one can perform a unitary transforma- 
tion which produces the following "shift" of the phonon operators: 



ak-^bk = ak + ai-^ bl. = al + 



The transformation 



is canonical: 



S = exp 



Vk 



Ok 



:i.2ii 



5"^ = exp 
and has the desired property: 

S~^akS = ak — 



V 



V, 



hid 



LO 



LO 



The transformed Hamiltonian is now: 

k k 

The eigenstates of the Hamiltonian contain an integer number of phonons (|?^k)) -The 
eigenenergies are evidently: 

|2 



Vk 



LO 



E = ^ rikhwLo - E 

k k 

This expression is divergent at it is often the case in field theory of point charges are consid- 
ered. A transformation of the type S has been of great interest in developing weak coupling 
theory as shown below. 



E. The "standard" theories 

1. Weak coupling via a perturbation theory 

For actual crystals a- values typically range from a — 0.02 (InSb) up to a ~ 3 to 4 (alkah 
hahdes, some oxides), see Table 1. A weak-couphng theory of the polaron was developed 
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originally by Frohlich [3|. He derived the first weak-coupling perturbation-theory results: 



Eo = -ahuLo (1.22) 

and 

m* = (1.23) 
1 — a/o 

Expressions (ll.22p and (11.231) are rigorous to order a. 

2. Weak coupling via a canonical transformation ["shift" -operators] 

Inspired by the work of Tomonaga on quantum electrodynamics (Q. E. D.), Lee, Low 



and Pines (LLP) 29j analyzed the properties of a weak-coupling polaron starting from a 
formulation based on canonical transformations (cp. the results of the subsection II D|) .As 
hown by them, the unitary transformation 

U = expl^^(p-^hkala^ -r^, (1.24) 

where P is a "c"-number representing the total system momentum allows to eliminate the 
electron co-ordinates from the system. Intuitively one might guess this transformation by 
writing the exact wave function in the form 

^totai H = exp (^^p ■ 1$) . 

It is plausible that the "Bloch" factor exp {i/hp ■ r) attaches the system to the electron as 
origin of the co-cordinates. After this transformation the Hamiltonian (I1.16p becomes: 

(p — hkalak) 

n = U-'HU = ^ ^ + J2 ^Loaj^ak + J^iV^a^ + V^al). (1.25) 

^"^'^ k k 

If, for the sake of simplicity, the case of total momentum equal to zero is considered, this 
expression becomes: 

" " + + ^ 4«k + $^(V^.ak + V:al). (1.26) 

k,k' k ^ ''^ k 

The first term of this Hamiltonian is the correlation energy term involving different values 
for k and k'.If one diagonalizes the second and the trird term of the Hamiltonian ( ll.26p 
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(this can be done exactly by means of the " shifted-oscillator canonical transformation" S 
fll.2ip ). the result of LLP is found. The expectation value of the first term is zero for the 
wave function ^lO) . Therefore one is sure to obtain a variational result. It is remarkable 
that merely extracting the k = k' term from the expression 

E /i^k ■ k'alal.,akakf 
k,k' " 

eliminates the divergency from the problem (cp. with the case nib oo) and is equivalent 
to the sophisticated theory by Lee, Low and Pines (LLP), which corresponds thus to neglect 
of the term fll.27p . The details of the LLP theory are given in Appendix 1. The explicit 
form for the energy is now 

^ = ~yi 1 — ^^^j r,2fc2 = -ahw. 

This self energy is no longer divergent. The divergence is elmininated by the quantum cut-off 
occurring at A; = y/2mhUJ-Lo/h. 

For the self energy the LLP result is equivalent to the perturbation result. The effective 
mass however is now given by 

m* =mh(l + ^^ , 

a result, which follows if one considers the case P ^0 and which is also exact for a — 0. 
However, the LLP effective mass is different from the perturbation result if a insreases. 

The LLP approximation has often been called "intermediate-coupling approximation". 
However its range of validity is the same as that of perturbation theory to order a. The 
significance of the LLP approximation consists of the flexibility of the canonical transfor- 
mations together with the fact that it puts the Frohlich result on a variational basis. 

To order a^, the analytical expressions for the coefficients are a 2; 21n(y2 + l)-|ln2- 



V2 



-0.0159 



2 

polaron mass 



L962 for the energy and | ln(v^ +l)-|ln2-^ + ^^ 0.02362763 for the 
3Q|. 



At present the following weak-coupling expansions are known: for the energy 3l|, l32 | 



° - -a - 0.0159196220a2 - 0.000806070048a^ - . . . , (1.28) 



LO 



and for the polaron mass 



30] 

— = ! + - + 0.02362763a2 + . . . (1.29) 
rrib 6 
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3. Strong coupling via a canonical transformation ["shift" -operators] 



Historically, the strong coupling limit was studied before all other treatments (Landau, 



Pekar 



21 



33|). Although it is only a formal case because the actual crystals seems to 



have a values smaller than 5, it is very interesting because it contains some indication of 
the intermediate coupling too: approach the excitations from the strong coupling limit and 
extrapolate to intermediate coupling is interesting because it is expected that some specific 
strong coupling properties "survive" at intermediate coupling. In what follows, a treatment, 
equivalent to that of Pekar, but in second quantization and written with as much analogy 
to the LLP treatment as possible is given. 

We start from the Frohlich Hamiltonian ( I1.16p . At strong coupling one makes the as- 
sumption (a "Produkt-Ansatz" ) for the polaron wave-function 



l<^>) = l^)l/) 



(1.30) 



where \(p) is the "electron-component" of the wave function ((v^lv^) = l).The "field- 
component" of the wave function |/) ((/|/) = 1) parametrically depends on \ip). The 
Produkt-Ansatz fll.30p — or Born-Oppenheimer approximation — implies that the electron 
adiabatically follows the motion of the atoms, while the field cannot follow the instantaneous 
motion of the electron. Frohlich showed that the approximation fll.30p leads to results, which 
are only valid for sufficiently large a — )■ oo, i. e. in the strong-coupling regime. A more sys- 
tematic analysis of strong-coupling polarons based on canonical transformations applied to 



the Hamiltonian fll.l6p was performed in Refs. 
The expectation value for the energy is now: 

,2 



34- 



W = M|-M + (/l 



I/) 



with 



We wish to minimize (H) , but also 



(/I ^/ia;Loa+ak + J](l4akPke^'^"- + \4*«U) I/) 

. k k 

has to be minimized. This expression will be minimized if |/) is the ground state wave func- 
tion of the "shifted" oscuUator-type Hamiltonian. As we can diagonalize this Hamiltonian 
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exactly: 



Ok 



LO 



2 I |2 



we can apply a canonical transformation similar to f ll.2ip : 



LO 



S = exp 



\ hi, It ^ 



flk 



which has the property: 



S-^ai,S = ak - 5-^4^ 
nujLo 



i Q - J 



VkPt 



The transformed Hamiltonian is now: 



(1.31) 
(1.32) 



;i.33) 



S ^ ^^WLoakOk + 5]^(VfcakPke'''"" + V^aJ^Pk) 

- k k 

Efc t \Pkf 
k k 

The phonon vacuum |0) provides a minimum: 



5 



{0\S 



-1 



/k^LoakCtk + J^(VfcakPke*^''" + ^aj^Pk) 



5|o) = -E 



LO 



Hence, the Hamiltonian ( ll.3ip is minimized by the ground state wave function 



S |0) = exp 
It gives the ground state energy 



^ \hWLO hWLO J 



|0). 



1.34) 



IT/ |2 I |2 

\Vk\ \pk\ 

hwi^o 



:i.35) 



which is still a functional of \ip). The functionals pk are different for differerent excitations. 
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a. Ground state of strong- coupling polarons For the ground state one considers a Gaus- 
sian wave function: 



\iP"is") = Cexp 



with a variational parameter Qq. 

/ \ 3 



-r 



h 



h 



3/4 



3/2 



dx exp 



rrib^o 2 

-X 



h 



3/2 



exp 



For the further use, we introduce a notation = (^^^jp) • Such a wave function is 
consistent with the locahzation of the electron, which we expect for large a. The kinetic 
energy in f ll.35p for this function is calculated using the representation of the operator 



p 



Is" I 



P 

2mh 



2mb 



C I d^rexp 



rribQo 2 

— r 

2h 



x(V^ + V^ + V^.)expf-^r^ 



2mb 
dy exp 



dxexp 
mb^lo 2 

-y 



2h 

rrib^l, 



2h 



2 

X 



V^exp 



TUb^O 



2h 



-X 



dz exp 



rrib^o 2 



—3 Cf I dxexp 



2mb 



2h 



-X 



TUb^O 



-xexp 



nib^l, 



2 

X 



2mb 

= 3^Ct , 
nib n 

6 — - — 



— oo 
oo 



2h 



dx { 1 ^"""^^ ) 



dx I 1 



— X I exp 

h 

nib^o 



2 ^ / nitUo 



h 



rrib^o 2 

-X 



h 
h 



-X 



,^A)^2 



rribQo 



-hQo. 
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The functional 



Pk"U 



h 



C"^ / d^r exp 



h 



r'^ + i 



h 



exp 



exp 



d^r exp 



-k • r- 



r + zk ■ r 



+ 



h 



r + 1 



2mb^c 



he 



The second term in (11 .351) is then 



-E 



IT/|2| |2 
l^fcl |Pk"ls"| 



V 



•lo 



(2 



d^k 



:7r 



fc2 Y y2mbUJi^o 



exp 



ah 



2ahhJi,o 



TX 



The variational energy fll.35p thus becomes 



3 hhJi^oa I VLq 



Putting 



one obtains 



a 



4 2^\l n 



9^0 



I ^0 _ 2a 

4^2 
^0 = t: — i^LO- 

9 TT 



fir 



4a2 



WlO 9 TT 



(1.36) 



(1.38) 



(1.39) 



Substituting (I1.39P in (ll.38p . we find the ground state energy of the polaron Eq (calculated 
with the energy of the uncoupled electron-phonon system as zero energy): 

/I 2\ «2 



3 4a2 hui^oa 2a , ^ - i " t 

4 9 vr y/n 3y/TT \3 3 / vr 

1 a;2 2^ 
Eq = -- — hwLo = -0.106a hWLO- 

3 TT 



(1.40) 



The strong-coupling mass of the polaron, resulting again from the approximation (ll.30p . is 
given [37| as: 

m* = 0.0200a^m6. (1.41) 
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More rigorous strong-coupling expansions for Eq and m* have been presented in the hterature 

Eq 



LO 



-0.108513a2 -2.836, (1.42) 



^ = 1 + 0.0227019a^ (1.43) 

Vflb 

The strong-couphng ground state energy (11.401) is lower than the LLP ground state energy 
for a > 10. 

b. The excited states of the polaron: SS, FC, RES In principle, excited states of the 
polaron exist at all coupling. In the general case, and for simplicity for P = 0, a continuum 
of states starts at hw-^o above the ground state of the polaron. This continuum physically 
corresponds to the scattering of free phonons on the polaron. Those "scattering. states" (SS) 



were studied in [391] anf for the first time more generally in [37|] are not the only excitations 
of the polaron. There are also internal excitation states corresponding to the excitations of 
the electron in the potential it created itself. By analogy with the excited states of colour 
centers, the following terminology is used. 

(i) The states where the electron is excited in the potential belonging to the ground state 
configuration of the lattice are called Franck-Condon (FC) states 

(ii) Excitations of the electron in which the lattice polarization is adapted to the electronic 
configuration of the excited electron (which itself then adapts its wave function to the new 
potential, etc. ...leading to a self-consistent final state), are called relaxed excited state 



(RES) |21|. 



c. Calculation of the lowest FC state The formalism used until now is well adapted to 
treat the polaron excitations at strong coupling. The field dependence of the wave function 
is (11.341) . For the FC state the pk are the same as for the ground state (ll.36p . Physically 
Pk tells us, to what electronic distribution the field is adapted. The electronic part of the 
excited wave function is 2p-like: 

|v."v') = C^"2p"^exp(^-^r2) (1.44) 



19 



with a parameter Qp, which is equal to : 



(V9"2p"|</5"2p") = C'"2p" / dzZ exp 



-z 



h 



1 TflijVlp 2 

dxexp I — X 



2p" \ rribflp 

h / 2 

3/2 



^2 



/I 



V 



3/2 



h 
h 



2 _ / "^b^p ^ 2mbQp 



We introduce still a notation 



3/2 



2 _ 2 / mbfip\ 
The FC state energy is, similarly to f ll.35p . 
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The kinetic energy term is 
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For the FC state, Qp = Qq. The second term in fll.45p is precisely fll.37p . 



— > r = ^V"i 



and the FC energy (11 .45^ becomes 
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The energy of the lowest FC state is, within the Produkt-Ansatz 40| : 



Efc = ^f^LO = 0.0354a^;ia;LO- 



[1A7) 



The fact that this energy is positive, is presumably due to the choice of a harmonic potential. 
The real potential the electron sees is anharmonic, and a bound state may be expected. 

d. Calculation of RES The electronic part of the excited wave function is (ll.44p with 
a variational parameter f2p, which is determined below. The variational RES energy is, 
similarly to ( ll.35p . 



E 



RES 



2p" 



2 iTr |2 I |2 

P I \ Pk"2p" 



2mb 



LO 



(1.48) 



Here the kinetic energy term is given by Eq. f ll.46p . The functional, which is now needed, is 
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Further, we substitute (ll.49p in the second term in the r.h.s. of Eq. fll.48p : 
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The variational energy ( I1.48P becomes 
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The energy of the RES is (see Refs. [37|, l39|): 

^RES = -0.042«2;^LO- 

The effective mass of the polaron in the RES is given [37(] as: 

uiIes = 0.621-^mb = 0.0200a^m6. 



;i.50) 



(1.51) 



The structure of the energy spectrum of the strong-couphng polaron is shown in Fig.[2l 
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FIG. 2: Structure of the energy spectrum of a polaron at strong coupling: — the ground 
state, -Eres — the (first) relaxed excited state; the Franck-Condon states {Eyq)- In fact, both the 
Franck-Condon states and the relaxed excited states lie in the continuum and, strictly speaking, 
are resonances. 

The significance of the strong-coupling large polaron theory is formal only: it allows to 
test "all-coupling" theories in the limit a — )■ oo. Remarkably, the effective electron-phonon 
coupling strength significantly increases in systems of low dimension and low dimensionality. 

Ji-. All-coupling theory. The Feynman path integral 



Feynman developed a superior all-coupling polaron theory using his path-integral formal- 



ism 



4l| . He studied first the self-energy Eq and the effective mass m* of polarons 41 1 
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Feynman got the idea to formulate the polaron problem into the Lagrangian form of 
quantum mechanics and then eliminate the field oscillators, ". . . in exact analogy to Q. E. D. 



. . . (resulting in) ... a sum over all trajectories . . . " . 
to the ground-state properties) is of the form (Ref. 



he resulting path integral (here limited 



(0,/3|0,0): 



T) exp 



13 



a 



-\r-a\ 



-drda 



r r 



(1.52) 
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where /3 = IjikBT). f ll.52p gives the amplitude that an electron found at a point in space 
at time zero will appear at the same point at the (imaginary) time /3. This path integral 
( I1.52P has a great intuitive appeal: it shows the polaron problem as an equivalent one- 
particle problem in which the interaction, non-local in time or "retarded", occurs between 
the electron and itself. Subsequently Feynman showed how the variational principle of 
quantum mechanics could be adapted to the path-integral formalism and he introduced a 
quadratic trial action (non-local in time) to simulate fll.52p . 

Applying the variational principle for path integrals then results in an upper bound for the 
polaron self-energy at all a, which at weak and strong coupling gives accurate expressions. 
Feynman obtained smooth interpolation between a weak and strong coupling (for the ground 
state energy). The weak-coupling expansions of Feynman for the ground-state energy and 
the effective mass of the polaron are: 

Eq 
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LO 
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In the strong-coupling limit Feynman found for the ground-state energy energy: 

Eq EsdM 



and for the polaron mass: 
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Over the years the Feynman model for the polaron has remained the most successful 
approach to this problem. The analysis of an exactly solvable ( "symmetrical" ) ID-polaron 



model 
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3, 



Monte Carlo schemes 



28 



43| and other numerical schemes 4J| demonstrate 



the remarkable accuracy of Feynman's path-integral approach to the polaron ground-state 
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energy. Experimentally more directly accessible properties of the polaron, such as its mobil- 
ity and optical absorption, have been investigated subsequently. Within the path-integral 
approach, Feynman et al. studied later the mobility of polarons 



45 



461 ]. Subsequently 



the path-integral approach to the polaron problem was generalized and developed to be- 
come a powerful tool to study optical absorption, magnetophonon resonance and cyclotron 
resonance 



47 



In Ref. 52|, a self-consistent treatment for the polaron problem at all a was presented. 



which is based on the Heisenberg equations of motion starting from a trial expression for 
the electron position. It was used to derive the effective mass and the optical properties of 
the polaron at arbitrary coupling. A variational justification of the approximation used in 



Ref. 



52| (through a Stiltjes continuous fraction) is reproduced in Appendix 2. 



5. On Monte Carlo calculations of the polaron free energy 



In Ref. 531], using a Monte Carlo calculation, the ground-state energy of a polaron was 



derived as Eq = lim^_!.oo AF, where AF = Fj^ — F2 with the free energy per polaron 



and FO = [3/ (2/3)] In (27r/3) the 
the actual computation in Ref. 



ree energy per electron. The value ptiuji^o = 25, used for 



531 ] ■ corresponds to T/Tq = 0.04 (Tp = hu^.n/kn; hui^o is 



the LO phonon energy). So, as pointed out in Ref. [5^], the authors of Ref. [5,3| actually 
calculated the free energy AF, rather than the polaron ground-state energy. 



To investigate the importance of temperature effects on AF, the authors of Ref. 



considered the polaron energy as obtained by Osaka 
polaron theory to nonzero temperatures: 



55| . who generalized the Feynman 
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This result is variational, with variational parameters v and w, and gives an upper bound 
to the exact polaron free energy. 
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The results of a numerical-variational calculation of Eq. f ll.57p are shown in Fig. |3l 
where the free energy —AF is plotted (in units of hwi^o) as a function of a for different 
values of the lattice temperature. As seen from Fig. [21 (i) — AF increases with increasing 
temperature and (ii) the effect of temperature on AF increases with increasing a. 
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FIG. 3: Contribution of the electron-phonon interaction to the free energy of the Feynman po- 
laron as a function of the electron-phonon couphng constant a for different values of the lattice 



temperature. Inset: temperature dependence of the free energy for a = 3. (From Ref. 



In Table IIT| the Monte Carlo results 53|, (AF)^^, are compared with the free energy of 
the Feynman polaron, (AF)p , calculated in 54|. The values for the free energy obtained 
from the Feynman polaron model are lower than the MC results for a < 2 and a > 4 (but 
lie within the 1% error of the Monte Carlo results). Since the Feynman result for the polaron 
free energy is an upper bound to the exact result, we conclude that for a < 2 and a > 4 the 
results of the Feynman model are closer to the exact result than the MC results of jssl . 
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TABLE II: Comparison between the free energy of the Feynman polaron theory, —(AF)p, and the 



Monte Carlo results of Ref. 



53l ]. — (AF)mC) for T/Td = 0.04. The relative difference is defined as 



A = 100 X [(AF)f - (AF)mc]/(AF)mc- (From Ref. [54]) 



a 


- (AF)f 


-(AF)mc 


A(%) 


0.5 


0.50860 


0.505 


0.71 


1.0 


1.02429 


1.020 


0.42 


1.5 


1.54776 


1.545 


0.18 


2.0 


2.07979 


2.080 


-0.010 


2.5 


2.62137 


2.627 


-0.21 


3.0 


3.17365 


3.184 


-0.32 


3.5 


3.73814 


3.747 


-0.24 


4.0 


4.31670 


4.314 


0.063 



6. On the contributions of the N-phonon states to the polaron ground state 



The analysis of an exactly solvable ("symmetric") ID-polaron model was performed in 



Refs. 
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, IsgI, [s^I . The model consists of an electron interacting with two oscillators possess- 
ing the opposite wave vectors: k and -k.The parity operator, which changes Ok and a_k (and 
also and a^j^), commutes with the Hamiltonian of the system. Hence, the polaron states 
are classified into the even and odd ones with the eigenvalues of the parity operator +1 and 
— 1, respectively. For the lowest even and odd states, the phonon distribution functions 
are plotted in Fig. 1, upper panel, at some values of the effective coupling constant A of the 
"symmetric" model. The value of the parameter 



{hky 



nibhuLo 

for these graphs was taken 1, while the total polaron momentum P = 0. In the weak- 
coupling case (A ~ 0.6) Wn is a decaying function of A^. When increasing A, Wn acquires 
a maximum, e.g. at = 8 for the lowest even state at A ~ 5.1. The phonon distribution 
function Wn has the same character for the lowest even and the lowest odd states at all 
values of the number of the virtual phonons in the ground state, (as distinct from the 



28 



higher states). This led to the conclusion that the lowest odd state is an internal excited 
state of the polaron. 



In Ref 



28 



the structure of the polaron cloud was investigated using the diagrammatic 
quantum Monte Carlo (DQMC) method. In particular, partial contributions of A^-phonon 
states to the polaron ground state were found as a function of for a few values of the 
coupling constant a, see Fig. 1, lower panel. It was shown to gradually evolve from the 
weak-coupling case (a = 1) into the strong-coupling regime {a = 17). Comparion of the 
lower panel to the upper panel in Fig. H] clearly shows that the evolution of the shape and 
the scale of the distribution of the iV-phonon states with increasing a as derived for a large 



polaron within DQMC method 



within the "symmetric" ID polaron model 



is in remarkable agreement with the results obtained 
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The "symmetric" polaron model (Devreese 1 964) 
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FIG. 4: Upper panel: The phonon distribution functions Wn in the "symmetric" polaron model 



5a], Fig. 23). Lower 



for various values of the effective coupling constant AatK = l,P = (from 
panel: Distribution of multiphonon states in the polaron cloud within DQMC method for various 
values of a (from [28,], Fig. 7). 
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F. Polaron mobility 



The mobility of large polarons was studied within various theoretical approaches (see Ref. 
for the detailed references). Frohlich 



59| pointed out the typical behavior of the large- 
polaron mobility 

/i cx exp{hui,of3), (1.59) 

which is characteristic for weak coupling. Here, /3 = l/fc^T, T is the temperature. Within 
the weak-coupling regime, the mobility of the polaron was then derived, e. g., using the 



60, 



6l| and starting from the LLP-transformation in Ref. 



Boltzmann equation in Refs. 

A nonperturbative analysis was embodied in the Feynman polaron theory, where the 
mobility fi of the polaron using the path-integral formalism was derived by Feynman et 
al. (usually referred to as FHIP) as a static limit starting from a frequency-dependent 
impedance function. For sufficiently low temperature T the mobility then takes the form 



45| 



w\ 3 3e 

V J AruhhulQaP 



;i.6o) 



where v and w are (variational) functions of a obtained from the Feynman polaron model. 

Using the Boltzmann equation for the Feynman polaron model, Kadanoff 63|| found the 
mobility, which for low temperatures can be represented as follows: 



3 e 
/i = — ( 



;i.6i) 



The weak-coupling perturbation expansion of the low-temperature polaron mobility as found 
using the Green's function technique |64] has confirmed that the mobility derived from 
the Boltzmann equation is exceedingly exact for weak coupling (a/6 ^ 1) and at low 
temperatures {ksT <^ ^lo)- As shown in Ref. (63|, the mobility of Eq. fll.60p differs by 
the factor of 3/(2/3/ia;Lo) from that derived using the polaron Boltzmann equation as given 
byEq. f lL6T|) . 



In the limit of weak electron-phonon coupling and low temperature, the FHIP po 



mobility of Eq. f ll.6Up differs b y th e factor of 3/(2/3^lo) from the previous result 



which, as pointed out in Ref. [45[ and in later publications (see, e.g., Refs. [58 



63 



60 



correct for /3 ^ 1. As follows from this comparison, the result of Ref. 45| is not valid when 



aron 



is 



T — )■ 0. As argued in Ref. 



45| and later confirmed, in particular, in Ref. 65| the above 
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discrepancy can be attributed to an interchange of two limits in calculating the impedance. 
In FHIP, for weak electro n-phonon coupling, one takes limn_j.o limQ,_j.o, whereas the correct 
order is limQ,_^o limn-^o is the frequency of the applied electric field). It turns out that 
for the correct result the mobility at low temperatures is predominantly limited by the 



absorption of phonons, while in the theory of 
the dominant contribution as T goes to zero 



'H 



65|. 



IP it is the emission of phonons which gives 



The analysis based on the Boltzmann equation takes into account the phonon emis- 
sion processes whenever the energy of the polaron is above the emission threshold. The 
independent-collision model, which underlies the Boltzmann-equation approach, however, 
fails in the "strong coupling regime" of the large polaron, when the thermal mean free path 



becomes less than the de Brog. 



be expected to be adequate 



45 



ie wavelength; in this case, the Boltzmann equation cannot 



66|. 



In fact, the expression f ll.60p for the polaron mobility was reported to adequately de- 
scribe the experimental data in several polar materials (see, e.g., Refs. 66l468| ) . Experi- 
mental work on alkali halides and silver halides indicates that the mobility obtained from 
Eq. f ll.60p describes the experimental results quite accurately 67|. Measurements of mo- 
bility as a function of temperature for photoexcited electrons in cubic ra-type Bii2Si02o are 
explained well in terms of large polarons within the Feynman approach 66|. The exper- 
imental findings on electron transport in crystalline Ti02 (rutile phase) probed by THz 
time-domain spectroscopy are quantitatively interpreted within the Feynman model GsJ. 
One of the reasons for the agreement between theory based on Eq. f ll.60p and experiment 
is that in the path-integral approximation to the polaron mobility, a Maxwellian distribu- 
tion for the electron velocities is assumed, when applying the adiabatic switching on of the 
Frohlich interaction. Although such a distribution is not inherent in the Frohlich interac- 
tion, its incorporation tends to favor agreement with experiment because other mechanisms 
(interaction with acoustic phonons etc.) cause a Gaussian distribution. 
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II. OPTICAL ABSORPTION. WEAK COUPLING 



A. Optical absorption at weak coupling [within the perturbation theory] 

At zero temperature and in the weak-coupUng Hmit, the optical absorption is due to the 
elementary polaron scattering process, schematically shown in Fig.O 




FIG. 5: Elementary polaron scattering process describing the absorption of an incoming photon 
and the generation of an outgoing phonon. 

In the weak-coupling limit {a ^ 1 ) th e polaron absorption coefficient was ffist obtained by 
V. Gurevich, I. Lang and Yu. Firsov 69|], who started from the Kubo formula. Their optical- 
absorption coefficient is equivalent to a particular case of the result of J. Tempere and J. T. 
Devreese (Ref. Isil), with the dyiiamic srtucture factor ^(q, u) corresponding to the Hartree- 



Fock approximation (see also 



70|, p. 585). At zero temperature, the absorption coefficient 



for absorption of light with frequency Q can be expressed in terms of elementary functions in 
two limiting cases: in the region of comparatively high polaron densities (^(f2 — u;lo)/C ^ 1) 



T(io) 



1 2V2Ar2/3a 



1 



and in the low-concentration region {h{Q — u)lo)/C ^ 1) 

1 2Ne^a {u - 1)^/^ 



e(a;-l). 



Tco) 



0(^-1), 



(2.1) 



(2.2) 



where u = Q/ui^o, co is the dielectric permittivity of the vacuum, n is the refractive index 
of the medium, N is the concentration of polarons and ( is the Fermi level for the electrons. 
A step function 

1 if w > 1, 
if a; < 1 
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reflects the fact tliat at zero temperature tlie absorption of light accompanied by the emission 
of a phonon can occur only if the energy of the incident photon is larger than that of a phonon 
(cj > 1). In the weak-coupling limit, according to Eqs. f l2.ip . fl2.2p . the absorption spectrum 
consists of a "one-phonon line" . 

At nonzero temperature, the absorption of a photon can be accompanied not only by 
emission, but also by absorption of one or more phonons. 



A simple derivation in Ref. 71 1 using a canonical transformation method gives the ab- 

n 

sorption coefficient of free polarons, which coincides with the result (12. 2p of Ref. |69| . 

B. Optical absorption at weak coupling [within the canonical-transformation 
method [71] (DHL)] 

The optical absorption of large polarons as a function of the frequency of the incident 
light is calculated using the canonical-transformation formalism by Devreese, Huybrechts 



and Lemmens (DHL) Ref. [TlJ. A simple calculation, which is developed below in full 
detail, gives a result for the absorption coefficient, which is exact to order a. 

We start from the Hamiltonian of the electron-phonon system interacting with light is 
written down using the vector potential of an electromagnetic field A if): 

= (p + -A it)) + 5^u;Loa^ak + + ^k«ke-"^") • (2.3) 
The electric field is related to the vector potential as 

Within the electric dipole interaction the electric field with frequency Q is 

E{t) = Ecos{nt) ^ (2.5) 

A = -^Esin(r]t) . (2.6) 
When expanding (p + ^A (t))^ in the Hamiltonian, we find 

2 „ ^2 
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where the first term is the kinetic energy of the electron, and the second term describes the 
interaction of the electron-phonon system with light 

Vt = A (t) ■ p = ■ p sin nt (2.8) 

Vt = Vsmnt, (2.9) 
^E-p. (2.10) 

Since A (t) does not depend on the electron coordinates, the term ^^^^2 (t) in (12. 7p does 
not play a role in our description of the optical absorption. The total Hamiltonian for the 
system of a continuum polaron interacting with light is thus 

Htot = Hpoi + Vt, 

where Hp^i is Frohlich's Hamiltonian fll.l6p . 

The absorption coefficient for absorption of light with frequency Q by free polarons is 
proportional to the probability P{^) that a photon is absorbed by these polarons in their 
ground state. 

Here is number of polarons, which are considered as independent from each other, Eq is 
the permittivity of vacuum, c is the velocity of light, n is the refractive index of the medium 
in which the polarons move, E is the modulus of the electric field vector of the incident 
photon. If the incident light can be treated as a perturbation, the transition probability 
P(fi) is given by the Golden Rule of Fermi: 

P{n) = 27iJ2^ {% \V\ f) (/ 1^1 $0) 6iEo + Q- Ej). (2.12) 

V is the amplitude of the time-dependent perturbation given by (I2.10p . The ground state 
wave function of a free polaron is |$o) and its energy is Eq. The wave functions of all possible 
final states are |/) and the corresponding energies are Ef. The possible final states are all 
the excited states of the polaron. The main idea of the present calculation is to avoid the 
explicit summation over the final polaron states, which are poorly known, by eliminating all 
the excited state wave functions |/) from the expression (I2.12p . 
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With this aim, the representation of the (5-function is used: 

fO 



1 r 

6{x) = — Re / dt exp [—i (x + ie) t] . 

J-oo 



This leads to: 



P{Q) = 2 Re 5^ r ($0 \V\ f) (/ \V\ %) exp + te + Eo - Ef)t] 

•' J — oo 

= 2Re^ / dtexp[-i{n + ie)t]{%\V\f){f\e'"'Ve~'"'\%). 



Using the fact that 
and the notation 

e'^*V^(0)e-'^* = V{t), 

dV{t) 
dt 

we find 



t [H, V{t)] 



P(fi) = 2Re/" dtexp[-i{n + ie)t\{%\V{0)V{t)\%) . (2.13) 

J — oo 

Defining 

R{n)= f dtexp[-i{n + ie)t]{%\V{0)V{t)\%) , (2.14) 

J —oo 

one has 

p(^]) = 2Rei^(^]). (2.15) 

Substituting fl210D to flCTj) . we find that 

R{n)=(^y f rfte-(^+-)*($o|(E-p(0))(E-p(t))|<|.o) (2.16) 

and hence 

P{Q) = 2(^) Re r ^ite-^(^+-)*(<|.o|(E-p(0))(E-p(t))|<|.o). (2.17) 

It is convenient to apply the first LLP transformation 5*1 (I4.187p . which eliminates the 
electron operators from the polaron Hamiltonian: 

H — y Ti = Si ^HpoiSi = 7io + 7ii : 
Hi = - — E k • k'alal^,ai,ai^,, 



k 



2mf, 
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where Ho can be diagonalized exactly and gives rise to the self-energy E = — Q;c<;LO,and 
Til contains the correlation effects between the phonons. The optical absorption will be 
calculated here for the total momentum of the system P = 0. 

In the LLP approximation the explicit form of the matrix element in fl2.16l) is 



($0 |(E • p (0)) (E • p $0) = (0 \S,'S^'E ■ pS.S^'E ■ p(t)^i^2| 0> 



(2.18) 



where 5*1 and S2 are the first fl4.187p and the second fl4.196p LLP transformations. The 
application of 5*1 gives: 

Using "H = S^^HSi, we arrive at S^^e^^^Si = e'^*. Further we recall S'f'^pS'i = P — 
^kaj^flk + p, where P = and p is set (see Appendix 1). This results in 



Then fl2.18p takes the form 
($o|(E-p(0))(E-p(t))|<|.o) = (o 



S^' ^ E ■ kala^ E • kaj,(t)ak(t) j S2 



Here the second LLP transformation is given by fl4.196p with 



fk 



v: 



V- 

(2.19) 
(2.20) 



and the vacuum is defined by ak|0) = 0. The calculation of the matrix element f l2.19p 
proceeds as follows: 







S2' E ■ kala^ E ■ kai(t)ak(t) j ^2 

^2"^ ^^E ■ kaj^flkj 8282^6''^' 828^^ kajcflkj 828^^e-''^'82 

82^ ^^E ■ kala^ ^2e*^2"'«S2t^-i • kaj^Okj 



) . (2.21) 



Further on, we calculate 



5*2 ^7182 = Hq + Hi 
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where 



Hq — 5*2 ^'HoS2 — 5*2 ^ 



^2 



Further we use 5*2 = ctk + /k: 



k 



l/k 



+ (c^LO + ^) (aL/k + ak/^) + [V, (ak + /k) + V,* (4 + Z^;) 

U2 



^ fcM t \Vk\' 



2mi, 



The last term can be calculated analytically: 



27r2 



1 \ 



Air f dk ( ) ' 

Jq \2mbUJLoJ 1 + 75- 



fc2 



2aui^o - 1 



TT 



2awLO 



1 + 



TT 



2mi,ajLO 

arctanx|^ = aui^o, 



The term 



will be neglected: 



Ho = -auj-Lo + X] (^^Lo + aL^k- 



i^i — 5*2 ^1^182 



Neglecting Hi, consistent with the LLP description, introduces no error in order a. There- 
fore ( I2.2ip becomes 







Ek E ■ k {ala^ + f^al + f*a^ + f^f*) e^^«* 
X Ek E ■ k (ala^ + f^al + f^a^ + e-^°* 



. 



(2.22) 
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For P = there is no privileged direction and E ■ k/^/^ = 0, fl2.22p reduces to: 



J]E-k/*ake^^°*J]E-k/k4e-^ 



Hot 



From the equation of motion for aj^ : 

da\.{t) 



dt 



Ho,al 



2mb 



it is easy now to calculate 



e " a^e 



« ujlo + 



2mh 



The matrix element f l2.19p now becomes 

($0 I (E ■ p (0)) (E ■ p m $o) = 5^ (E ■ k)^ /^/k exp 



2mh 



The transition probability fl2.13p is then given by the expression 



P{n) = 2 Re ■ k)' /k/k /° dt exp 



i ( i7 + - wlo 



2m;, 



2mb 
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Using fl2.20p . we obtain 



27re^ 



[E ■ k)^ 



y 

2m^ J 



X 



(2 



1 



fc2 

2mi, 



,6 n 



V \2mbUJ-Lo 



2mt 



2n 



-1 



dk 



2mh 
1 



2mbULo 



) 



1 



X 



-6 



2m,5ajLO / 



A;2 



8eW 



r5 



2mbUJi,oJ 



X 



e 



- 1 



n 



- 1 



where 







n 

UJhO 



1 if 
if 



n 
n 



> 1 
< 1 



finally takes 



(2.23) 



The absorption coefficient (12. lip for absorption by free polarons for a 
the form 

The behaviuor of rp(f2) f l2.23p as a function of Q is as follows. For Q < wlo there is no 
absorption. The threshold for absorption is at = WLO-From Q = cjlo up to = I^lo, Tp 
increases to a maximum and for Q > |a;Lo the absorption coefficient decreases slowly with 
increasing Q. 

Experimentally, this one-phonon line has been observed for free polarons in the infrared 
absorption spectra of CdO-films, see Fig. [61 In CdO, which is a weakly polar material with 
a ~ 0.74, the polaron absorption band is observed in the spectral region between 6 and 20 
fim (above the LO phonon frequency). The difference between theory and experiment in the 
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wavelength region where polaron absorption dominates the spectrum is due to many-polaron 
effects. 




1 d , I , I , I , 

10 20 30 40 50 

X (|jm) 



FIG. 6: Optical absorption spectrum of a CdO-film with the carrier concentration n = 5.9 x 10^^ 
cm~^ at T = 300 K. The experimental data (solid dots) of Ref. [t^] are compared to different 
theoretical results: with (solid curve) and without (dashed line) the one-polaron contribution of 

n n 

Ref. [69] and for many polarons (dash-dotted curve) of Ref. |51| . 



III. OPTICAL ABSORPTION. STRONG COUPLING 

The absorption of light by free large polarons was treated in Ref. [47[ using the polaron 
states obtained wihtin the adiabatic strong-coupling approximation, which was considered 
above in subsection IIE3 



It was argued in Ref. 47j, that for sufficiently large a (a > 3), the (first) RES of a polaron 
is a relatively stable state, which can participate in optical absorption transitions. This idea 
was necessary to understand the polaron optical absorption spectrum in the strong-coupling 
regime. The following scenario of a transition, which leads to a "zero-phonon" peak in the 
absorption by a strong-coupling polaron, can then be suggested. If the frequency of the 



incoming photon is equal to 



ite = ; = 0.065« wlo, 
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then the electron jumps from the ground state (which, at large coupling, is well-characterized 
by "^''-symmetry for the electron) to an excited state ("2p"), while the lattice polarization 
in the final state is adapted to the "2p" electronic state of the polaron. In Ref. {4^! 
considering the decay of the RES with emission of one real phonon it is demonstrated, that 
the "zero-phonon" peak can be described using the Wigner-Weisskopf formula valid when 
the linewidth of that peak is much smaller than ^lo- 
For photon energies larger than 

^RES + ^LO, 

a transition of the polaron towards the first scattering state, belonging to the RES, becomes 
possible. The final state of the optical absorption process then consists of a polaron in 
its lowest RES plus a free phonon. A "one-phonon sideband" then appears in the polaron 
absorption spectrum. This process is called one-phonon sideband absorption. 

The one-, two-, ... K-, ... phonon sidebands of the zero-phonon peak give rise to a 
broad structure in the absorption spectrum. It turns that the first moment of the phonon 
sidebands corresponds to the PC frequency Qpc- 

_ Efc - Eq 2 
^FC = z = 0.141q; cjlo- 

h 

To summarize, the polaron optical absorption spectrum at strong coupling is characterized 
by the following features (T = 0): 

a) An intense absorption peak ("zero-phonon line") appears, which corresponds to a 
transition from the ground state to the first RES at ^res- 

b) For Q > Qres + i^lo, a phonon sideband structure arises. This sideband structure 
peaks around fipc- 



The qualitative behaviour predicted in Ref. j47|, namely, an intense zero-phonon (RES) 
line with a broader sideband at the high-frequency side, was confirmed after an all-coupling 



expression for the polaron optical absorption coefficient at a = 5, 6, 7 had been studied [48|. 

In what precedes, the low-frequency end of the polaron absorption spectrum was dis- 
cussed; at higher frequencies, transitions to higher RES and their scattering states can 
appear. The two-phonon sidebands in the optical absorption of free polarons in the strong- 
coupling limit were numerically studied in Ref. 



73|. 
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The study of the optical absorption of polarons at large coupling is mainly of formal 
interest because all reported coupling constants of polar semiconductors and ionic crystals 
are smaller than 5 (see Table 1). 



IV. ARBITRARY COUPLING 

A. Impedance function of large polarons: An alternative derivation of FHIP [tJ ] 

a. Definitions We derive here the linear response of the Frohlich polaron, described by 
the Hamiltonian 

H=^ + Y1 ^ka,^«k + ^^(Vfcake^'^-'- + l^^aj.e"*^-'-), (4.1) 
k k 

to a spatially uniform, time-varying electric field 

En(t) = ^oexp(zfit)e,. (4.2) 

This field induces a current in the x-direction 

Mt) = ^En(t). (4.3) 

The complex function Z{Q) is called the impedance function. The frequency-dependent 
mobility is defined by 



For nonzero frequencies (in the case of po 
one defines the absorption coefficient 48| 



arons the frequencies of interest are in the infrared) 



r(fi) = ^Re-^, (4.5) 

where eo is the dielectric constant of the vacuum, n the refractive index of the crystal, 
and c the velocity of light. In the following the amplitude of the electric field Eq is taken 
sufficiently small so that linear-response theory can be applied. 

The impedance function can be expressed via a frequency-dependent conductivity of a 
single polaron in a unit volume 
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using the standard Kuho formula (cf. Eq. (3.8.8) from Ref. 70|): 

.2 1 



+ 



e^^"(b.(t),jj)rft. 



(4.7) 



In order to introduce a convenient representation of the impedance function, we give in 
the next subsection a definition and discuss properties of a scalar product of two operators 
[cf. y, Chapter 5]. 

b. Definition and properties of the scalar product For two operators A and B (i.e., 
elements of the Hilbert space of operators) the scalar product is defined as 



{A,B) = / rfA((e^''Mt)5). 
Jo 



(4.8) 



The notation (e^^'^A''') is used in order to indicate that the operator e'^^^'" acts on the operator 
. The time evolution of the operator A is determined by the Liouville operator L: 



with a commutator [H, A] , wherefrom 

A{t) = e'^'A{0) = e'^'^^A{0)e-'^'/\ 
The expectation value in ( 14. 8 p is taken over the Gibbs' ensemble: 

{A) = Tt{po{H)A) 
with the equilibtium density matrix when the electric field is absent 



(4.9) 



(4.10) 



(4.11) 



Po{H) = e-^^/Tr(e 



(4.12) 



One can show that (14. 8 p defines a positive definite scalar product with the following prop- 
erties 



(i) 



n 



(iii) 



{A,B) 
{A, LB) 
{A, LB) 



and [cf. Eq. (5.11) of [75|] 



(B^A^), 

{LA,B), 
1 

h 



(iv) {A,Br = {B,A) 



(4.13) 
(4.14) 
(4.15) 

(4.16) 
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(4.17) 



Demonstration of the property (14.131) . Starting from tlie definition (14. 8p and using 
flirTOj) . we obtain 

Jo 

Substituting here ( 14. lip with f l4.12p . one finds 

{A,B) = fdXTi [e-(^-")^Ate-^^i?] /Tr(e-^^). 
Jo 

Change of the variable A' = /3 — A allows us to represnt this integral as 



iAB) 



dX'Ti 



/Tr(e 



Further, a cyclic permutation of the operators under the trace Tr sign gives 

{A,B) = fdXTi [e-(^-^)^i?e-^^At] /Tr(e-^^) 
Jo 

= r^A (e^^5e-"^At) = f d\ {{e^'^'^ B) X^) 
Jo Jo 



According to the definition (14. 8p . this finalizes the demonstration of f l4.13p . 



Demonstration of the property ( 14.141) . Starting from fl4.17p and using (14. 9p . we 
obtain 

{A,LB)= / dX{e^^A^e~^^LB) 
Jo 

= i l^dX (e^^Ate-^^ {HB - BH)) . 
^Jo 

A cyclic permutation of the operators under the average (•) sign gives 



1 

(A, LB) = - dX (e^^A^-^^HB - He^^A^-^^B) . 
^Jo 

Using the commutation of H and e^^^ , one finds 

1 

{A, LB) = - dX {e^^A^He^^^B - e^"HA^e-^^B) 
= i l^dX (e^^ {A^H - HA^) e-"^5> 
= (e^^ (HA - AH^ g-^^E^ 

^/^.A(e-(i[^,^])^-i.) 



(4.18) 

(4.19) 
(4.20) 
(4.21) 

(4.22) 
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B 



With the definition f l4.9p . this gives 

{A, LB) = j^X (^e^^ {LA)^ e^^^fi^ = . j^X ^''^ {LA)^ 
According to the definition fl4.8p . this finahzes the demonstration of fl4.14p . 

Demonstration of the property ( I4.15P . Starting from fl4.20p and performing a cychc 
permutation of the operators under the average (•) , we find 

1 

{A, LB) = ^ dX (e^^ {A^H - HA^) 

Further we notice that 



dX 

consequently, 

^^^'^^—nj/\ dx ^ 
'h \^Jq dX 

= _^ (|e^^Ate-^^i?|^) = liiA^B - e^"" A^e-^"" B)) 

= ^Tr [e-^^ {A^B - e^"" A^-^"" B)] /Tr(e-^^) 

= ^Tr le-f^'A^B - Ah-^^^B] /Tr(e-^^). (4.23) 
n 

Further, a cychc permutation of the operators in the second term under the trace Tr sign 
gives 

{A, LB) = -^Tr [e-^"A^B - e'^^BA^] /Tr(e-^^) 
= -^Tr [e-^^ {A^B - BA^)] /Tr(e-^^) 
= l{A^B-BA^) = l{[A\B]). 



Thus, the property (14.150 has been demonstrated. 
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Demonstration of the property (14.161) . We start from the represntation of the 
scalar product fl4.17p and take a complex conjugate: 

io Jo 
A cyclic permutation of the operators under the average (•) sign gives then 

{A,By= fdX{e'^"B^e-'^"A)= TdA (e^^^te'^^A) 
Jo Jo 

dX{{e^''^B^) A) . (4.24) 



According to the definition fl4.8p . this finalizes the demonstration of fl4.16p . 

The above scalar product allows one to represent different djTiamical quantities in a 
rather simple way. For example, let us consider a scalar product 



A, 



z-L 



-B 



(4.25) 



r-oo rP 



B 



(4.26) 



POO pp 

-I / dte''' / d\ (e^'^^A^e-'^'B) 
Jo Jo 

POO pP 

-i j dte''^ / d\ {e^^A^e'^^e-'^^/'^Be'"^/'') 
Jo Jo 

POO pf} 

-I dte"' t/A (e^^*/''+^^Ate-^^*/''-^^5) 
Jo Jo 

POO 

-I / dte'^' / dX (^e'm-im|n^^^-^H(t-ixn)/n^^ 
Jo Jo 

POO pP 

I dte''' d\{e'^^'-'^^^A^B) 
Jo 

dX{A\t-thX)B{0)) . 



oo i-p 

izt 



—I I dte 

'0 JO 



(4.27) 



(4.28) 



c. Representation of the impedance function in terms of the relaxation function The 
impedance function is related to the relaxation function 

1 



$ (z) = (z) 



-X 



Z-L 

where x is the velocity operator, by the following expression: 

1 



(4.29) 



(4.30) 



[z = n + ie,e> 0). 
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Demonstration of the representation ( 14.291) . Apply fl4.27p to the relaxation 
function entering (14.301) : 



dX 



and perform the integration by parts using the formula 



izt 



e^'f{t)dt=—f{t)e 

z 

= r/(0) + - 



t=0 

oo 



+ - 

Z 



df it) 



dm ^^.t 

dt 



dt 



e'^'dt 



izt ^^^i{t—ihX)L^ 



x)x)dt = -{ {e^^^x) x) + 



oo a 
izt ^ 



dt 



^i(t—ihX)L 



c) x) dt 



- { (e'^^^x) x) - - / e*^* {L (e^(*-^''^)^i;) x) dt. 



This allows us to represent the relaxation function in the form 



^izt ^^^it-^HX)L.^^ 



/3 



1 i 
— + - 



dX { [e^^^x] x) + - I dX I e''' {L (e^(*-'^^)^x) x) dt 

/3 



dX 



^izt ^^iit-ihX)L^~^ 



rribZ z Jq 

where the expression fl4.55p is inserted in the first term. Further on, the integral over A is 
taken as follows: 



/3 



dX {L (e'^'-'^^^^x) x) 



L 



j{t^inx)L^^^ 1 X 



^^iLt ^^hpL^ -x)x) 



h 
1 

h 

\{{e''^''x {t)-x{t))x) 



1 



Tr [e-^^ (e^^x (t) e"''^ - x (t)) x] 
Tr [(x (t) e-^^-e-^^x (t)) x] 



^;^r^Tr[e-^^(xx(t)-x(t))x] 
'xx (t) — X (t) x) 



(x(t) ,x) 
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Hence, 



oo 

izt 



When setting z = Q + ie with e — > +0, we have 



$ (z) = ^ / dte''' {[x (t) , x]) . (4.31) 



oo 



For f2 7^ 0, we can set 

$ + = — - / e^(^+'^)* ([x (t) , x]) dt. (4.33) 



Multiplying $ (f2 + ie) by ie"^, we find that 

ze'* (fi + ^e) =ie^(—l--^ H e^(^+^-)* ([i (t) , i]) rft 
yTTih \l nil J fj 

= 2 — - + |- / e'(^+^^)* ( [i; (t) , x] ) c^t 



e2 1 /.oo 



mf,f2 h^l Jq 
where the electric current density is 

Jx €X. 

Substituting further (14.340 in (I4.30p . we arrive at 

1 e2 , 1 



^ / e^(^+-)*([j.(t),jj)rft, (4.34) 



1 /"°° 



Z(i7) m^i^ 

what coincides with the expression of the impedance function fl4.6p through a frequency- 
dependent conductivity given by the Kubo formula fl4.7p . q.e.d. 

d. Application of the projection operator technique Using the Mori-Zwanzig projection 
operator technique (cf. |75|], Chapter 5), the relaxation function fl4.29p 

^ (z) = { X, -X 



z-L 

can be represented in a form, which is especially convenient for the application in the theory 
of the optical absorption of polarons. 

The projection operator P {Q = 1 — P) is defined as 

PA = (4.35) 
X 

with A an operator and 

X=(x,x). (4.36) 
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The projection operator Q = 1 — P projects an operator onto the space orthogonal to the 
space containing x. Here we give some examples of the action of the projection operators: 



Px = x, Qx = (1 - P)x = 0; 

x{x,x) x{iLx,x) ix{x,Lx) 



Px 



X 



IX 



X 



X X 
Qx = (1 — P)x = X. 

i (x, Ok) x{iLx,ak) ix{a]f_,Lx) ix 

X 



{[x,x]) = 0, 



Pa 



X X 
Qttk = (1 - P)ak = Ok 



X 



The projection operators P and Q are idempotent: 



P'A 



x(x,A) \ 

2 a x{x,A) 



PA: 



(4.37) 
(4.38) 
(4.39) 
(4.40) 
(4.41) 



X X 
Q"^ = {l- Pf = l-2P + P"^ = 1- P = Q. 

The Liouville operator can be identiacUy represented as L = LP+LQ. Then the operator 
^3^in the relaxation function (14.291) can be represented as follows: 



z-L z-LQ-LP 
We use the algebraic operator identity: 

1 _ 1 1 1 

X + y X x^ X + y 



with X = z — LQ and y = —LP : 

1 



z — L z — LQ z — LQ z — L 
Consequently, the relaxation function (I4.29P takes the form 

1 



^iz) 



X \ + \ X 



' LP- ' 



z — LQ J \ z — LQ z — L 
The first term in the r.h.s. of fl4.42p simplifies as follows: 
1 



-X 



X, 



z-LQ 



-X 



- + 1-lq + \lqlq + ... 



z z^ 







(■ l-\ 


x^ 













because Qx = 0. Using the quantity fl4.36p we obtain: 

1 



z-LQ 
50 



X 



(4.42) 



(4.43) 



The second term in the r.h.s. of (14.421) contains the operator 

P^—x 
z — L 

which according to the definition of the projection operator P f l4.35p can be transformed as 

1 X I 1 \ X 

P X = - (x, X = -$ (z) . (4.44) 

z-L x\z-Lj X 

It is remarkable that this term is exactly expressed in terms of the sought relaxation function 
f lCT]) . Substituting f 05]) and f iOIj) in f 02]) . we find 



2;$ (z) = X + ( i 



z \ ' z-LQ X 



^ L- 1 $ (z) =^ 



' :L-]<l>(z) 



z-LQ X 

. z — LQ + LQ ^ x \ ^ , , 

^ + 1 ^' in ' "^^-^^ 

z-LQ X. 



= X + 1 a;, 
2;$ {z)=x + 



1 + 



i, Lx) 1 
X X\z-LQ 



L-\^(z) 



Lx 



X 

LQ 



<^(z) 



Introducring the quantity 



O 



(x, Lx) 



X 



and the function called the memory function 



^{z) = - (x,LQ \-pzLx ) , 

X V z-LQ 



(4.45) 



(4.46) 



(4.47) 



we represent (14.451) in the form of the equation 

[z-0-^iz)]'^{z)=x. 

A solition of this equation gives the relaxation function $(-2) represented within the Mori- 
Zwanzig projection operator technique: 



X 



(4.48) 



z-0-^{z) 

The memory function (I4.47P can be still transformed to another useful form. First of all, we 
apply the property of a scalar product (I4.14p : 
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j:{z) = -(lx,q^-— 

X \ z- LQ 



Lx 



-((P + Q)Lx,Q^-— 
X\ z- LQ 



Lx 



-(PLX,Q 
X\ z-LQ 



Let] + - (qLx,Q \— 

J X\ z-LQ 



Lx 



(4.49) 
(4.50) 
(4.51) 



For any two operators A and B 

{PA,QB) 



f X {x, A) 


x{x,By 


V X ' 


X . 



{x,A) 
X 

X 



x,B)- 
[(x,S)-( 



{x, x) {x, B) 
X 

x.B)] =0, 



therefore the first term on the r.h.s. in fl4.49p vanishes, and we obtain 

S(z) = - (qLx,Q \r—Lx 

X \ z-LQ 



(4.52) 



In this expression, the operator Q can be represented in the following form, using the 
fact that Q is the idempotent operator: 

1 



Q 



z-LQ 



Q 



- + ^LQ + \lQLQ 



z z^ 



-Q + \qLQ + \qLQLQ + ... 

Z Z^ Z"^ 

^-Q + ^QLQ"" + \qLQ''LQ^ + ... 



z 



z^ 



- + \qlq + \qlqqlq 



Q- 



1 



z z^ 
1 



Q 



(4.53) 



z — LQ z — QLQ 

A new Liouville operator can be defined, C = QLQ, which describes the time evolution 
in the Hilbert space of operators, which is orthogonal complement of x. Substituting then 
(I4.53P with the operator C into (14.521) . we bring it to the form, which will be used in what 
follows. 

S(2) = - (qLx, -^QLx] . (4.54) 
For the Hamiltonian (11.161) we obtain the following quantities: 



X = [x,x) 



Px .J 

— , iLx 
rrib 



TJibh 



{Px, Lx) 



52 



Using fl4.15p . we find 



and 



= ^^ l = m,-{[x,x]) = 0. 



Substituting (14.551) and (14.561) in (I4.48p . one obtains 

1 1 



The operator 



Lx = L 



1 



rrib z — 



Px _ 

1 



rribh 



Px 



t *k-r\ 



rUb 



Lx 



n, -i— ^ 



nib 



(4.55) 
(4.56) 

(4.57) 



does not depend on the velocities. Therefore, multiplying both parts of (jl? 
taking into account (14.391) and (14.411) . we obtain 



(4.58) 
with Q and 



QLx = -— V UVka^e""^-^ - V^ale-'^^-n 

TUb 



what allows us to represent the memory function in the form 



Yj{z) = - ( QLx, — ^—tQLx 




=—yyKKVkv:, ^ ' • 4.59 

^ k k' \ I^zlOk'e +aj^,e 

In transition to (I4.59P we have used the property of the amplitude (I1.17p : = —Vk and 
taken into account that according to the definition (14.80 . the first operator enters a scalar 
product in the hermitian conjugate form. Introducing the operators 

&k = ake^''^6i = aie-*'^^ 
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we represent the memory function as 



E{z) = —Y.Yl kXVkV:> f (6k + ^(&k' + bl)) . (4.60) 



k k 

„ikr\ „ ^ik-r 



We notice that Q6k = Q(«ke ) = = 6k- It will be represented through the four 

relaxation functions: 

S(^) = —Y.Y. kXVkV:, [$+ + (^) + $-^7(z) + $+k7(^) + $- + (^)] , (4.61) 



" k k' 



(4.62) 

Kyi^) = (^k,^^6k'), (4.63) 



2; 


-C 




1 


z 


-C 




1 


z 


-£ 




1 


z 


-£ 



'^kkt(^) = ^k,-^6LV (4-65) 



There exist relations between the above relaxation functions. For example, the relaxation 
function f l4.63p . takes the form 

'^kk^(^) = (&k,T^^k' 







z-C 

dte'^' (e^^*6k(0),6k'(0)). 



Then the complex conjugate of this relaxation function: 

(^)]* = ' / dt^"'*' (e^''*6k(0),6k'(0))* 

POO 

= t dte-^^**(6k'(0),e^^*6k(0)), 
Jo 

where the property (I4.16p has been used. The property (I4.13p gives 

[$-j^7(^)]* = t / dte-''*' (6k'(0), e^^*/''6k(0)e-'^*/^) 
Jo 



00 



tj rfte-^^**(e*^*/''6^(0)e*^*/^6i,(0) 
[dte~-'' (6L(t),6U0))=-$:k' 
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wherefrom it follows that 

%yi^) = -[Kyi-^nY- (4-66) 

Similarly, the relation 

\-.tW = -K.T(-.*)]* (4.67) 

is proven. 

e. Memory function In this subsection we indicate which approximations must be 
made in the calculation of the relaxation functions in order to obtain the FHIP results 
for the impedance function. Consider the relaxation function fl4.62p : 



^ dte'''(^e'^%liO),bliO)'^ =-tJ^ dte"' (bi{t),bl{0) 



where b]^{t) = e*^*6|^(0), and perform a partial integration: 



\ (&ic(0),&U0)) + -J^i^'''' {^Cbi{t)M^)) (4.68) 



z 

Here we supposed that 



^ 'fe|c(0),&U0)) - rdte-' to(t),&L(0) 



hm e-* 6L(t),6L,(0) =0. 

t — >oo V / 

In the second term in fl4.68p . 

£6[(t),6U0)) = Ue^''%,bi) = (QLQe'Q'^Q%,bl{0) 
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because Qbl = bl and = Q. Further on, we have 



cbi{t),bi{o) 



So, we find from fl4.68p 



1 



*kk' W = (&lc(0),&L(0)) - z I dte^^^Lbi{t),bl{0) 



6k(t),&U0) 



The first term in the r.h.s. of this expression can be represented as follows: 
I {bU0),bU0)) = \[dx{{e'''b,) bi) = \ Q^X bl 



hL 



1. (e-i..e-.. - i.. 



\lblXb^-lb^bl 



zh 



L L 



— i^j dte'HyJbl- j dtble'^-'bi, 
6k(t),&U0) 



I 

zh 



dt 



Substituting it in the r.h.s. of f l4.69p . we find 



zh 

i 



dt 



zh 



dte 



izt 



j^l dt{i-e-^'){ &k(t),&Uo) 
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In a similar way one obtains 



• poo 

W = ^j^ dt{l- e-*) ([6k(t),&k'(0)]) . (4.71) 

Inserting the relaxation functions (14 .70 p . (I4.7ip . (I4.66P and (I4.67p . we find the memory 
function (I4.6ip 



1 i r°° 

" k k' -"^ 



\ ( 






&k(t),&Uo)" 



+ ([&k(t),&k'(0)]) 
^ k k^ ^/iJo 



Im f(k(t),&U0)l) + ([&k(t),&k'(0)]) 



wherefrom 



1 f°° 

S(z) = -/ rft(l-e*"*)lmF(t) 



with 



m = -AeE^-^^^^^^' {(f&k(t),6Uo)l) + ([&k(t),6k'(o)]) 



nibh 



(4.72) 



(4.73) 



k k' 



/. Derivation of the memory function To calculate the expectation values in Eq. (14. 73 p . 
we shall make the following approximations (cf. Ref. 76|). The Liouville operator C, which 
determines the time evolution of the operator ^^(t) = e*^*6|^(0), is replaced by Lph + Lp, 
where Lph is the Liouville operator for free phonons and L^t^ is the Liouville operator for the 
Feynman polaron model 4l| . The Frohlich Hamiltonian appearing in the statistical average 
(•) is imilarly replaced by Hph + if i;', with Hph the Hamiltonian of free phonons and Hp the 
Hamiltonian of the Feynman polaron model. With this approximation, e.g., the average 



K{t)bl{0)) = (ak(t)aUo)) (e^'^-^We-^'^'-') = 4,^' («k(t)at(0)) (e^'^'^We"^'^-') . (4.74) 
The time evolution of the free-phonon annihilation operator f l4.10p . 

akit) = e^^p'^^/'^akC"*^'''^*/'^ = exp {iu]^a^akt) akcxp [-iuj]ia^ai^t) 
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accorting to (14. 9p is 

.dak(t) 

-I' — 
dt 



= Wkexp (i^kOkC^k^) [oil^^a-k] o-kexp (^—iui^a^akt) 
= — WkGxp (iWkCtk'^k^) Ok exp (— iaJkOjI^akt) =^ 
= exp {—icoift) Ok. 



Similarly, 



al(t) = exp (iukt) a 



k- 



Hence, we have 



a]i{t)alj = exp (-iwk^) {^Cikolj = exp {-iukt) (l + ala^j = exp {-iLo^t) [1 + ^(wk)] , 

where n(cjk) = [exp(/3^(X'k) — 1] ""^ is the average number of phonons with energy hw]^. 

The calculation of the Fourier component of the electron density-density correlation func- 
tion <||e«k r(t)g-jk r^^ f|4.74l) for an electron described by the Feynman polaron model is 



given below following the approach of Ref. 76] . 
We calculate the correlation function 



/ ik.r(t) -ik.r(r)\ _ ^ ' ^ ^ ) (a 

\ I Tr(e-^^^) ^ ^ 

with the Feynman trial Hamiltonian 

ifF = — + -^ + -x(r-r.)'. (4.76) 

Here, r (t) denotes the operator in the Heisenberg representation 

r (t) = e^^^^re"^^^. (4.77) 

We show that the correlation function {^e«k r(t)g-«k r(r)^^ depends on (r — t) rather than on 
t and r independently: 

^gik.r(t)g-ik-r(r)^ ^ _V 



Tr(e-^^^) 
Tr(e-/5^^) 
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where a = t — t. 

The normal coordinates are the center-of-mass vector R and the vector of the relative 
motion p: 



R 



mT+mfTf 



m+rrif 

p = r — Tf 



The inverse transformation is: 



Ri f 

m+inf' 

J m+rrif' 

The same transformation as (14.791) takes place for velocities: 



P 




m+mf 

From f l4.79p we derive the transformation for moments 

p _ P _|_ '"Z Pp 
— — — \ ~r — ; m m f 

m m+inr m+mf J 

P/ _ P m Pp 
— ; — \ mm ^ 

m J m-\--m f m-\-m j 



' m+m J 

mt y-, 

Pf = — r—^ — Po 

^ J m+m f 

The Hamiltonian fl4.76p then takes the form 



with the masses 



m + mf 



(4.79) 



(4.80) 



(4.81) 



(4.82) 



p2 p2 1 - 
Hf = — + ^ + -p^y (4.83) 



TflTfi f 

M = m + mf,p= (4.84) 



^ = J-. (4.85) 



and with the frequency 

The Cartesian coordinates and moments corresponding to the relative motion can be in the 
standard way expressed in terms of the second quantization operators: 

P.i = -«(^)'"(Q-Cj). (4.86) 
(i = 1.2,3) 
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In these notations, the Hamiltonian f l4.83p takes the form 

HF = ^ + i2^(cjC, + l). (4.87) 

j=l V / 

Using fl4.87p . we find the operators in the Heisenberg representation, (i) for the center-of 
mass coordinates 

= + ~ ^'^^^^ + ^'^^^^ " ^^^^^ 



-^T^Pj {-2ih) = Xj + 
2Mh ^ ^ ^ ^ M 



X, {cj) = X, + (4.88) 

(ii) for the operators Cj and Cj 

Cj (a) = Cje-'^", C] [a) = Cje'^". (4.89) 
Using the first formula of fl4.80p we find 

r (a) = R (a) + -^p (a) 

1/2 



+ (Ce-- + CV-). (4.90) 



We denote 



1/2 / f,^2 n1/2 



M \2fiQ) ~ \2M^fin) 



/ hm f \ 

^ J ' m+mt / 



m+m f 

a 



r (a) = R + — P + a (Ce"'^'' + Ch'^^j . (4.91) 



Therefore, we obtain 

e'^ '' = exp (zk ■ R) exp [iak ■ C + iak ■ C"^) 

3 

= JJ^exp (ikjXj) exp {iakjCj — iakjCj 
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cr 



-zk ■ R - i—k ■ P - iak- ( Ce"^^" + C"^e^^" 
M 



a 



exp ( -zk ■ R - i—k ■ P ) exp ( -lak ■ Ce"^^" - zak ■ C"^e*^" 



JJ^ exp y-ikjXj - '''^^jPjj 6xp y-iakjCjC '^'^ - iakjCje'^" 



The disentangling of the exponents is performed using the formula 

gA+B ^ e^Texp (^j^ dXe-^^Be^^^ . 
In the case when [A, B] commutes with both A and B, this formula is reduced to 



We perform the necessary commutations: 



cr 



1 . o cr 



i fcj , i — kj Pj — -^k-— [Xj , Pj 



1 1 — 1- — ^cr, 



^ 2/2 

-a /c^. 



Cj, Cj 

^ 2(2 

-a k^ 



^ 2)2 

-a kj 



^ 2(2 

--a k, 
2 ^ 



c«k-r _ ^ik-R^iak-Ct iak-C^-ia2fc2 
c — c c c c ^ , 



g-ik.r(a) ^ ^-ik-Rg-i|^k.Pg-iak.Cte>"-g-iak.Ce-'"-gi||^a-|a2fc2 



ik-r -ik-r((T) _ -i-§k-P Mk-Ct iak-C -jak-Cte'""^ -jak-Ce-'"'' j||^(T-a2fc2 

It follows from Eq. (14.941) that when [A, B] commutes with both A and 5, 



Using (14.951) . we find 



iaki-Ci -iaki-CU^^" -iak^-CU^^'^ iak^-Ci kaki-C.-iakrCU^^"] 

o J J p J J — p 3 p J J pY J J ' ■'J J 



-iakj -Cje^""^ iakj -Cj a'^k'^e'^" 
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Herefrom, we find 

The correlation function then is 



since the variables of the center-of mass motion and of the relative motion are averaged 
independently. 

^gmk-Ct(l-e'"'^)giak-C(l-e-^^'^)\ _ /giQ-CtgiQ*-C^ 



^g-/3ftnCt.CgiQ CtgjQ*.C^ 

~ Tr (e"-^^^^^'^) 

with 

Q = a - e^^'^) k. 
Let us consider the auxiliary expectation value 

AQ*C\ — \ J 



1 A (^g)'^ A w 



V / n=0 m=0 



oo 

m=0 
oo 

m=0 

where |m) are the eigenstates of (C^C). The operators C act on these states as follows: 

C \m) — \/rn\m — 1) , 
C|0) = 0. 

Therefore, we find 

(m I (C^)" C"| m) = m (m — 1) . . . (m — n + 1) = '■ — -r for n < m, 

\ i\ / I / [m — ny. 

(m I (C^)'' C"| m> = for n > m. 
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oo 



-f)Klm 



{m — n)V 



and 



oo (_i)n|g|2n oo 



ml 



n=0 ^ ' m=n 
°° I i\".|^|2n oo 



(m — n)! 



n=0 m=n 

n=0 

oo 

= E 



-/ShQm 



m 



n 



n=0 

oo , „ I ,2„, 



fc=0 



k + r 
\ n 



1 



n=0 



1 y. (-ingr" 

n=0 ' \ ) 



1 



n+l 



exp 



In particular, ior Q — Q* — we have 

V /I - e-^'^" 

As a result, the expectation value /e^^C'^e*'^**" ) is 



e'^^^e^«*^) = exp[-n(f7)|g|^]. 



with 



Using this result, we obtain the expression 



giQCtgiQ*.C 



exp [-n (Q) Q • Q* 



= exp 
= exp 



-4n(Q) a2A;^sin^ -fia 
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The expectation value (e ^li^''^'^ is 



/dPexp(-/j^-^^k-P 

/dPexp(-/3^) 
/ dP exp 



2Af/3 2M/3 



/ciPexp(-/3^) 



Collecting all factors in Eq. ( I4.97P together, we find 



{ 



gik-rg-ik-r(o-) 



) = exp 



An (n) a^k"^ sin^ ( -no] -a^P (l-e 
2M 2M/3 ^ ^ V2 / V 



3ik-r(t) -ik-r{(T+i) 



> 



with the function 



1 - exp(ifit) + 4n (fi) sin' 



where 



M= — mfe,l2 = t;wLO,a =- — . 



(4.100) 



(4.101) 



(4.102) 



According to fliTTSj) . 



,ik-r(t)^-ik-r(cr+t)^ _ ^^ik-r^-ik-r(o-) 



e e 



Taking a = — t in f l4.100p we finally find the Fourier component of the electron density- 
density correlation function ^e«k-r(t)g-?.k-r^ which enters Eq. f l4.74p : 



^gik.r{t)g-ik.r^ = exp [-PD(-t) 
Finally, the correlation functions in f l4.73p reduce to 

\{t)bl{0)^ = [1 + n{uk)] exp (-zwkt) exp [-k^D{-t)] , 

(&k(t)&k(0)) = 0, 
(&k(0)fek(t)) = 0. 

Inserting these equations into Eqs. fl4.72p and (14. 73 p . one obtains 

S(z) = - / dt(l-e''')lmS(t) 
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(4.103) 



(4.104) 



with 

T cu^ 2 ^ 2itai2t J -[l + ^(wk)]exp(-iwkt)exp[-A;2D(-t)] 

lmS{t) = -> K\Vk\ lm< 

^bh ^ y +n(wk) exp (-^Wkt) exp [-A;2D(t)] 

Using the property D{—t)* =^ D(t) for real vaues of t,one obtains 

— Imexp (— icuk^) exp [— t)] } = Im exp (iwkt) exp [— A;^D(— t)*] } 

= Imexp (iwk^) exp [—k'^D(t)'j } 



and consequently 
2 

= — ^X^^^ iHil^exp [-k'^D{t)] {[1 + n(c<;k)] exp (icukt) + n(c<;k) exp (-zcuk^)} 

k 



(4.105) 



Owing to the rotational invariance of \Vk\ and a;k,we can substitute in f l4.105p 



The resulting expression for 
2 



S{t) 



k'^ l^fcl^exp [—k'^D{t)'j {[1 + n{coy^)] exp (iwkt) + n(co'k) exp (— ic^kt)} 

k 

(4.106) 

is identical with Eq. (35) of FHIP j45|. In the case of Frohlich polarons, taking into account 
f lTTTl) . Eq. (|4A06|) simplifies to 



"^6 / 3v27r I jf-n{uLo) exp (-iwLO^) 



(4.107) 



B. Calculation of the memory function (Devreese et. al. 48l|) 



Upon substituion of f l4.30p and f l4.57p into Eq. (14. 5p . we find the absorption coefficient 



Re 



neoc 



1 



1 e' 



^]-s*(^]) 



Im 



_f]-S(f]) 



Im 



neoc mfe [[n- Re E(n)]^ + [Im S(fi)]^ 
1 ImS*(^]) 
ne^^[r^-ReS(fi)]^ + [ImS(fi)]^ 
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r(fi) 



1 



ImS(r]) 



(4.108) 



This general expression was the strating point for a derivation of the theoretical optical 



absorption spectrum o 
Devreese et al. in Ref. 



a single large polaron, at all electron-phonon coupling strengths by 
. The memory function E(f2) as given by Eq. (14. 104p with fl4.107p 



48 



contains the dynamics of the polaron and depends on a, temperature and f2. Following the 
notation, introduced in Ref. |45[, 



E(^]) 



n 



we reresent Eq. ( I4.108P in the form used in Ref. 



48 



1 



(4.109) 



(4.110) 



According to fl4.104p and fl4.109p . 

lmx{^)=lm dtsm(nt)S{t), Rex(^^)=Im/ dt [1 - cos{nt)] S (t) . (4.111) 
Jo Jo 

In the present Notes we limit our attention to the case T = (/3 — ?■ oo). It was demonstrated 



in Ref. 



48 



that the exact zero-temperature limit arises if the 



imit /3 — )■ oo is taken directly 



in the expressions (14.1 lip (see Appendices A and B of Ref. 48|). As follows from f l4.107p . 

3/2 



S{t) = 

Accorting to (I4.10ip 



a 



3727? 



exp [iui^ot) [D{t)] 



-3/2 



(/3^oo) 



(4.112) 



. ht 

''2M 



+ a^[l- exp(zfit)] (/3 ^ oo). 



Using the Feynman units (where h = IjCJlo = 1 and mi, = 1), we obtain from (14.1021) : 

,,2 _ „„2 

2" 



M 



f\2- 2 IV — 



W 



fl = v,a 



and consequently 

^ , , 1 — 
= 2—' 

with 



:i-e™*)-i 



2\vJ 



1 /W\ 2 



2 \v 



{i?(l - e'^*) - zt} (/3^oo) (4.113) 



R 
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and according to (14.1121) 



S{t) 



2a fv\'^ 



3i/7r \w 

From (14.1 lip one obtains immediately 



e'' [R{1 - e™*) 



it] 



-3/2 



(/3^oo) 



(4.114) 



Im 
Im 



dt- 



sin(f2t)e* 



[i?(l 



□rat 



)-zt] 



3/2- 



[1 - cos(fit)] e'^ 



(4.115) 
(4.116) 



45| . However, to study 



[i?(l-e™*) -it]^/^' 

In the limit /3 — )■ oo the function Imx(i7) was calculated by FHIP 
the optical absorption to the same approximation as FHIP's treatment of the impedance, 
we have also to calculate Rex(^) and use this result in (I4.110p . The calculation of Rex(^), 
which is a Kramers-Kronig-type transform of Imx(f2), is a key ingredient in Ref. J48|. The 
details of those calculations are presented in the Appendices A, B and C to Ref. [48 1. 

Developing the denominator of both integrals on the right-hand side of (I4.115P and 
(I4.116p . the calculations are reduced to the evaluation of a sum of integrals of the type 



Im 



°° sin(fit)e*(^+"'^)* 

3/2+n • 



Im 



{R - it) 



dt 



cos(^]^)e*(^+'^'')* 



{R - it) 



3/2+n 



(4.117) 



In Appendix B to Ref. 48| it is shown how such integrals are evaluated using a recurrence 
formula. For Imx(i7) a very convenient result was found in 48|: 



Imx(fi) 



2a / f \ 3 



3 

X \Vt 



w 



3/2 1 



-1)' 



n=0 

nv 



n+1/2 ^-\n~l-nv\ 



(2n + l)...3-l 

^1 + sgn(i7 - 1 



nv 



(4.118) 



This expression is a finite sum and not and infinite series. FHIP gave the first two terms of 
( HTT8|1 explicitly. 

Using the same recurrence relation it is seen the analytical expression (see Appendix B to 
Ref. |[48|]), which was found for Re x(^) is far more complicated. To circumwent the difficulty 
with the numerical treatment of Rex(^), the corresponding integrals in (I4.117P have been 
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transformed in 



48] to integrals with rapildy convergent integrands: 
[1 - cos((]t)] e*(i+"'^)* 



Im 



dt- 



{R - itf^'' 



r(n + |; 



dx 



X In 



'l + nv + x) 



fi^ — (1 + nf + x) 



1/2 



(4.119) 



The integral on the rig ht-hand side of f l4.119p is adequate for computer calculations. In 



Appendix C to Ref. 



48 



some supplementary details of the computation of fl4.119p are given. 
Another analytical representation for the memory function (I4.104p was derived in Ref. \(\ - 



C. Discussion of optical absorption of polarons at arbitrary coupling 

At weak coupling, the optical absorption spectrum fl4.108p of the polaron is determined 
by the absorption of radiation energy, which is reemitted in the form of LO phonons. For 
a > 5.9, the polaron can undergo transitions toward a relatively stable RES (see Fig. [7]). 
The RES peak in the optical absorption spectrum also has a phonon sideband-structure, 
whose average transition frequency can be related to a FC-type transition. Furthermore, 
at zero temperature, the optical absorption spectrum of one polaron exhibits also a zero- 
frequency "central peak" [~ 5(r2)]. For non-zero temperature, this "central peak" smears 
out and gives rise to an "anomalous" Drude-type low-frequency component of the optical 
absorption spectrum. 



For example, in Fig. [7] from Ref. [48|], the main peak of the polaron optical absorption 
for a = 5 at n = 3.51u;lo is interpreted as due to transitions to a RES. A "shoulder" at 
the low-frequency side of the main peak is attributed to one-phonon transitions to polaron- 
"scattering states". The broad structure centered at about Vt = G.Scjlo is interpreted as a 
FC band. As seen from Fig. U\ when increasing the electron-phonon coupling constant to 



a=6, the RES peak at = 4.3wlo stabilizes. It is in Ref. [48| that the all-coupling optical 
absorption spectrum of a Frohlich polaron, together with the role of RES-states, FC-states 
and scattering states, was first presented. 

Recent interesting numerical calculations of the optical conductivity for the Frohlich 
polaron performed within the diagrammatic Quantum Monte Carlo method {t?], see Fig. |8l 
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10 



FIG. 7: Optical absorption spectrum of a polaron calculated by Devreese et al. 481] a = 5 and 6. 
The RES peak is very intense compared with the FC peak. The frequency Vt/uJi^o = ?; is indicated 
by the dashed lines.) 



fully confirm tlie essential analytical results derived by Devreese et al. in Ref. |48l for a ^ 3. 
In the intermediate coupling regime 3 < a < 6, the low-energy behavior and the position of 
the RES-peak in the optical conductivity spectrum of Ref. [77| follow closely the prediction 
of Ref. 48|. There are some minor qualitative differences between the two approaches in the 
intermediate coupling regime: in Ref. jl^l, the dominant ("RES") peak is less intense in the 
Monte-Carlo numerical simulations and the second ("FC") peak develops less prominently. 
There are the following qualitative differences between the two approaches in the strong 
coupling regime: in Ref. |77| . the dominant peak broadens and the second peak does not 
develop, giving instead rise to a flat shoulder in the optical conductivity spectrum at a = 6. 
This behavior has been tentatively attributed to the optical processes with participation of 
two jrsl or more phonons. The above differences can arise also due to the fact that, within 
the Feynman polaron model, one-phonon processes are assigned more oscillator strength and 
the RES tends to be more stable as compared to the Monte-Carlo result. The nature of the 
excited states of a polaron needs further study. An independent numerical simulation might 
be called for. 
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FIG. 8: Left-hand panel: Monte Carlo optical conductivity spectra of one polaron for the weak- 
coupling regime (open circles) compared to the second-order perturbation theory (dotted lines) 



for a = 0.01 and a = 1 and to the analytical DSG calculations 



48[ (solid lines). Right-hand 



panel: Monte Carlo optical conductivity spectra for the intermediate coupling regime (open circles) 



compared to the analytical DSG approach 
phonon thresholds. (From Ref. [t^]-) 



48l | (solid lines). Arrows point to the two- and three- 
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In Fig. IHl Monte-Carlo optical conductivity spectrum of one polaron for a = 1 compares 
well with that obtained in Ref. jrS] within the canonical-transformation formalism taking 
into account correlation in processes involving two LO phonons. The difference between the 
results of these two approaches becomes less pronounced when decreasing the value of a = 1 
and might be indicative of a possible precision loss, which requires an independent check. 
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FIG. 9: One-polaron optical conductivity Iiea{uj) for a = 1 calculated within the Monte Carlo 



781] (solid 



approach [77(] (open circles) and derived using the expansion in powers of a up to 
curve) . 

The coupling constant a of the known ionic crystals is too small {a < 5) to allow for 
the experimental detection of sharp RES peaks, and the resonance condition Q = i?eS(f2) 



48|. Nevertheless, for 3 < a < 5.9 the 



cannot be satisfied for a < 5.9 as shown in Ref, 
development of RES is already reflected in a broad optical absorption peak. Such a peak, 
predicted in Ref. was identified, e. g., in the optical absorption of Pr2Ni04.22 in Ref. 



791]. Also, the resonance condition can be fulfilled if an external magnetic field is applied; 
the magnetic field stabilizes the RES, which then can be detected in a cyclotron resonance 
peak. 
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1. Sum rules for the optical conductivity spectra 



In this section, we analyze the sum rules for the optical conductivity spectra obtained 
within the DSG approach 481 with those obtained within the diaErrammatic Monte Carlo 
calculation [77[. The values of the polaron effective mass for the Monte Carlo approach are 
taken from Ref. [28 1. In Tables 3 and 4, we represent the polaron ground-state Eq and the 
following parameters calculated using the optical conductivity spectra: 



Mo 
Ml 



Re a iuj) du. 



uRea (co) du. 



(4.120) 
(4.121) 



where cjmax is the upper value of the frequency available from Ref. 



rzl, 



Mo 



TC 



2m* 



Rea (u) du. 



where m* is the polaron mass, the optical conductivity is calculated in units 



(4.122) 



, m IS 



measured in units of the band mass rrih, and the frequency is measured in units of cjlo- The 
values of Wmax are: Uma.x = 10 for a = 0.01, 1 and 3, Wmax = 12 for a = 4.5, 5.25 and 6, 
oJma^ = 18 for a = 6.5, 7 and 8. 
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Table 3. Polciron pEirameters obtained within the diagrammatic Monte Ccirlo 

approach 



a 




^*(MC) 








0.01 


0.00249 


1.0017 


1.5706 


0.634 


-0.010 


1 


0.24179 


1.1865 


1.5657 


0.65789 


-1.013 


.3 


0.67713 


1.8167 


1.5280 


0.73123 


-3.18 


4.5 


0.97540 


2.8742 


1.5219 


0.862 


-4.97 


5.25 


1.0904 


3.8148 


1.5022 


0.90181 


-5.68 


6 


1.1994 


5.3708 


1.4919 


0.98248 


-6.79 


6.5 


1.30 


6.4989 


1.5417 


1.1356 


-7.44 


7 


1.3558 


9.7158 


1.5175 


1.2163 


-8.31 


8 


1.4195 


19.991 


1.4981 


1.3774 


-9.85 
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Table 4. Polaron parameters obtained within the path-integral approach 



a 


^(DSG) 

U 


(Feynman) 


^(DSG) 
u 


± / 


^(Feynman) 



0.01 


0.00248 


1.0017 


1.5706 


0.633 


-0.010 


1 


0.24318 


1.1957 


1.5569 


0.65468 


-1.0130 


3 


0.69696 


1.8912 


1.5275 


0.71572 


-3.1333 


4.5 


1.0162 


3.1202 


1.5196 


0.83184 


-4.8394 


5.25 


1.1504 


4.3969 


1.5077 


0.88595 


-5.7482 


6 


1.2608 


6.8367 


1.4906 


0.95384 


-6.7108 


6.5 


1.3657 


9.7449 


1.5269 


1.1192 


-7.3920 


7 


1.4278 


14.395 


1.5369 


1.2170 


-8.1127 


8 


1.4741 


31.569 


1.5239 


1.4340 


-9.6953 



The parameters corresponding to the Monte Carlo calculation are obtained using the 
numerical data kindly provided by A. Mishchenko. The comparison of the zero frequency 



moments M^^^ and M)^^^^ with each other and with the value 7r/2 corresponding to the 



r(DSG) 



sum rule 



80| 



shows that 



(MC) 



(DSG) 



h / Re cr [ijS] dcu = — 

2m* ^ ^ 2 

is smaller than each of the differences ? — M, 



(4.123) 



(MC) 



' 2 



(DSG) 







which appear due to a finite interval of the integration in (14. 12011 . (14.1211) . 
We analyze also the fulfilment of the ground-state theorem jsil 



Eo («) - ^0 (0) 



/ Lo ilea {uj,a ) auj 

TT Jq a' 



(4.124) 



using the first-frequency moments M^^'^^^ and m}^^*^\ The results of this comparison are 
presented in Fig. [TOl The dots indicate the polaron ground-state energy calculated using the 
Feynman variational principle. The solid curve is the value of Eq (a) calculated numerically 
using the optical conductivity spectra and the ground-state theorem with the DSG optical 
conductivity [48| for a polaron, 



E, 



(DSG) 



a) 



rdo^ r^Re^(DSG)^^^^/)^^_ (4_^25) 
Jo «' Jo 
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The dashed and the dot-dashed curves are the values obtained usine; M[^^^\a) and 



m}^*"-* (a), respectively: 



a' 



TC 



da 



(DSG) 



(MC) 



a) 



3 rda^ r^\Rea^^''\uj,a')duj = -- r da 

Jq a' Jo TT Jq 



a 

(MC) 



a' 



a' 



(4.126) 
(4.127) 



As seen from the figure, E^^^^'^^ (a) to a high degree of accuracy coincides with the vari- 



ational polaron ground-state energy. Both eI^^^^ (a) and E, 
due to the integration over a finite interval of frequencies. However, E^^^^^ (a) and E'q^'"^ (a) 
are very close to each other. Herefrom, a conclusion follows that for integrals over the finite 
frequency region characteristic for the polaron optical absorption (i. e., except the "tails"). 



n{MC) 



a) differ from E^^'^^ {a) 



the function Eq^'^^ (a) fl4.127p reproduces very well the function E^ 



(DSG) 



a] 



o 




FIG. 10: The ground-state theorem for a polaron using different data for the optical conductivity 



n 

spectra, DSG from Ref. [48] and MC from Ref. 



771 ]. The notations are explained in the text. 
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D. Scaling relations 



1. Derivation of the scaling relations 

The form of the Frohhch Hamiltonian in n dimensions is the same as in 3D, 

2 

H=^ + Y. ^^kaj,ak + (Vkake*^-^ + V^aj.e"^'^ '-) , (4.128) 
^ k k 

except that now all vectors are n-dimensional. In this subsection, dispersionless longitudinal 
phonons are considered, i.e., = wlo, and units are chosen such that h = nib = i^lo = 1- 
The electron-phonon interaction is a representation in second quantization of the electron 
interaction with the lattice polarization, which in 3D is essentially the Coulomb potential 
1/^- l^kl^ is proportional to the Fourier transform of this potential, and as a consequence 
we have in n dimensions 

\V.f = J^, (4.129) 

where L„ is the volume of the n-dimensional crystal. Note that \V^^f, where k is an {n — 1)- 
dimensional vector, can be obtained from |Vk|^, where k = ^k, /c^^ is an n-dimensional 
vector, by summing out one of the dimensions explicitly: 



\Vi.\' = Y\^^\'- (4.130) 



Inserting Eq. (14.1291) into Eq. f l4.130p . we have 



Replacing the sum in Eq. (14.1311) by an integral, i.e. 



L 



n-l 



Y^hl dkn, (4.132) 



we obtain 



27r (j, , ,,\("-l)/2 



(P + kl 



(4.133) 



Since 



(zf + x^Y r (p) 
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(4.134) 



we have 



1 r(|)r(^) _ y?r(n^) 



where F (x) is eh T function. Inserting Eq. fl4.135p into Eq. (14. 1331) . we obtain 

"^'^ P (— ) 

In 3D the interaction coefficient is well known, |Vk|^ = 2\f2'nal L^^k^ , so that 



(4.135) 



(4.136) 



As = 2^2: 



vra. 



(4.137) 



Inserting Eq. f l4.137p into Eq. fl4.136p . we immediately obtain 

r(l) 



A2 = 2V2' 



Tca- 



V2^ 



Tca. 



Applying Eq. fl4.136p n — 2 times, we further obtain for n > 3 

n-l 



(4.138) 



(2 A)""' n r (I) 



i=2 



n-2 

nr(i) 



-A, 



r(i) 



(4.139) 



So, the interaction coefficient in n dimensions becomes 



2n-3/2^{n-l)/2r (n_l) ^ 



(4.140) 



Following the Feynman approach j41|, the upper bound for the polaron ground-state 
energy can be written down as 



E = Eo-\ua-{S-So),. 

/3— s>oo p 



(4.141) 



where S is the exact action functional of the polaron problem, while Sq is the trial action 
functional, which corresponds to a model system where an electron is coupled by an elastic 
force to a fictitious particle (i.e., the model system describes a harmonic oscillator). Eq is 
the ground-state energy of the above model system, and 

/ Fe^opr (t) 



(4.142) 
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As indicated above, the Frohlich Hamiltonian in n dimensions is the same as in 3D, 
except that now all vectors are n-dimensional [and the coupling coefficient |Vk|^ is modified 
in accordance with Eq. fl4.140p ]. Similarly, the only difference of the model system in 
n dimensions from the model system in 3D is that now one deals with an n-dimensional 
harmonic oscillator. So, directly following |4l|], one can represent lim^^oo {S — So)q //3 as 

lim l(5-5o)o = A + (4-143) 

/3-s>oo p 

where 

POO 

A = Y^ \Vy.f / (e*''-['-W-^(°)]>ge-*dt, (4.144) 



B = ^-^ -j^ ([r(t)-r(0)]0^e-*dt, (4.145) 

w and V are variational parameters, which should be determined by minimizing E of Eq. 
f l4.14ip . Since the averaging (...)q in Eq. fl4.144p is performed with the trial action, which cor- 
responds to a harmonic oscillator, components of the electron coordinates, Tj {j = 1, ...,n). 



in (^e^Hrit)-rm^^ separate |4lj: 



^ik.[r(t)-r(0)]\^^ ^ J-j- ^gifcj[r-j(i)-r,(0)]^^ _ (4.146) 



For the average (^e*'^^^'"^*^*^ ^^^^^^)o^ Feynman obtained 41 1 



'0 

where 



^gifc,[r,(t)-r,(0)]^^ = ^-k^Doit)^ (4.147) 



W = is' + ^(1 (4.148) 
Inserting Eq. 04.1461) with Eq. fl4.147p into Eq. 04.1441) . we obtain 

POD 

A= e'*rft Virkl'e-'^'^^W. (4.149) 
Jo k 

Inserting expression O4.140p for |Vk|^ into Eq. 04.149P and replacing the sum over k by an 
integral [see f l4.132p ]. we have 



n 
2; 



2 J Jo J A;"-i {27vy 

\ roD p POO Q—k^l 



2n-3/2^(n-l)/2r I ^ T f T 6'^ j^O W /.n-l^^fc ^ ^^^^^^^ 
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where dQn is the elemental solid angle in n dimensions. Since the integrand in Eq. fl4.150p 
depends only on the modulus k of k, one have simply J dQn = Qn with 



So, we obtain for A the result 

2-1/2^-1/2^ (n_l) ^ 



A 



r(t) 

2-3/2^ (n_l) 



r(t)- 



(4.151) 



2-1/2^-1/2^ (n_l) ^ 



r(t) 



-.dt 

2v//^o(t) 

(4.152) 



Like in Ref. 



4l|, B can be easily calculated by noticing that 



;[r(t)-r(0)f>„ = 5^([r,(t)-r,(0)f>„ = 5^ 

i=i 



i=i 

n 



,ik-[r(t)-r(0)p 



k=0 



SO that 



B 



nw 


y 


-m;2) 




2 




nw 


y 


-ti;2) 




2 




nw 


y 


-^i;2) 




2 




n {v 


2 _ 


^2) 



J]2Do {t) = 2nDo (t) 



Do (t) e-'"'dt 

2^2 



(4.153) 







f 2 — w"^ 



2t;3 



dt 



} -1 V V 

W V V — w 

+ 



2f 2 f 2 



1 



2t>^ \w f + w 



(4.154) 



Inserting Eq. f l4.140p with A and 5, given by Eqs. f l4.152p and f l4.154p . together with 
the ground-state energy of the model system 4l| (an isotropic n-dimensional harmonic 
oscillator), 

Eo 



n{y — w) 



(4.155) 



into Eq. (14.1410 . we obtain 
E 



niv — w) n ff 2 — w"^] 



2-3/2p |n^^^ ^ 



n[v ~ w] 



4f 

r(^)« 



r(f) 



oo 



dt 



4.V 2V2r (f ) 



oo ^-t 



Do (t) 



-.dt. 



(4.156) 
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In order to make easier a comparison of E for n dimensions with the Feynman result 41 1 
for 3D, 

3 (f — w)'^ 



EsD (a) 



1 

--=a 
JO 



dt. 



it is convenient to rewrite Eq. f l4.156p in the form 



EnD (a) 



n 
3 



'?>{v-wf 1 3y?r(^^) 



a 



-.dt 



(4.157) 



(4.158) 



2nr(f) Jo 

It is worth recalling that the parameters w and v must be determined by minimizing E. 
Thus, in the case of Eq. (14.1581) one has to minimize the expression in the square brackets. 
The only difference of this expression from the r.h.s. of Eq. (14.1571) is that a is multiplied 

by the factor 



3v^r (^) 



(4.159) 



2nr(f) 

This means that the minimizing parameters w and v in nD at a given a will be exactly the 
same as those calculated in 3D for the Frohlich constant as large as a„a: 



VnD (a) = f3D {dnO) , WnD («) = W3D {cLnOt) 



(4.160) 



Therefore, comparing Eq. (14.1581) to Eq. (14.1571) . we obtain the scaling relation 821-184 1 



EnD («) 



(4.161) 



where is given by Eq. (14.1591) . As discussed in Ref. [82|], the above scaling relation is 
not an exact relation. It is valid for the Feynman polaron energy and also for the ground- 
state energy to order a. The next-order term (i.e., a^) no longer satisfies Eq. (I4.16ip . The 
reason is that in the exact calculation (to order a^) the electron motion in the different 
space directions is coupled by the electron-phonon interaction. No such a coupling appears 
in the Feynman polaron model [see, e.g., Eq. (I4.146P ]: and this is the underlying reason for 



the validity of t 
In Refs. 82 



le scaling relation for the Feynman approximation. 



8J], scaling relations are obtained also for the impedance function, the 



effective mass and the mobility of a polaron. The inverse of the impedance function Z (cu) 
is given by 

^ ' (4.162) 



Z (w) w — S {oj) 
where the memory function S (w) can be expressed as [85] 

S(z) = - / dt{l-e''')\mS{t) 
z Jo 
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(4.163) 



with z = uj + iO"'" and 

S(t) = Y^ 2kl \Vv\^ e-^'^^T (t) , (4.164) 



k 



T{t) = [l + n (1)] e* + n (1) e"**, (4.165) 



I vf 

1 _ e"™* + 4n (t;) sin^ 



Here, /3 is the inverse temperature and n (u) is the occupation number of phonons with 
frequency u (recall that in our units cjlo = !)• 

As implied from Eqs. f l4.162p and f l4.163p . scaling of S (u) and Z [u) is determined by 
scaling of S (t). For an isotropic crystal, since \Vkf, D {t) and T {t) do not depend on the 
direction of k, one can write ^? I'^kl^ e-^^'^^T (t) = fc| |Vk|^ e-'^'^^'^T (t) = ... = 
Ek K l^kl^ e-^'^WT [t) , so that 

S{t) = -yk'' l^kl' e-'^'^WT (t) . (4.167) 

k 

Inserting expression fl4.140p for |\4;|^ and replacing the sum over k by an integral, we have 

?^ (27r) y Jo 

n (27r)" 1(1)70 ^' 

= \[^^-^^^ r k'e~''^^'^T it) dk. (4.168) 
\ IT nV (f ) Jo 

In particular, for 3D one has from Eq. fl4.168p 

(a; t) = ^« /t^e-'^'^WT (t) rffc. (4.169) 
Stt Jo 

For nD, Eq. (14.1681) can be rewritten is the form 

^.D (a; t) = ^ 'f^iy c. /"^ A;^e-^-Wr (t) c^fc 
Stt 2r2r (I) Jo 



2^2 
-a 



POO 

.„« / k^e-^'''^'^T{t)dk. (4.170) 



So, the only difference of the expression for S^b (t) from 6*30 (t) is that a is multiplied by 
a„. Since for the minimizing parameters w and v, which enter D (t), scaling is determined 
by the same product a with a„ [see Eq. fl4.160p ]. we can write 

SnDia;t) = S3r,{ana;t) , (4.171) 
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so that [84 1 



(4.172) 



and 
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ZnD (a; = (a„a; u) . (4.173) 
The polaron mass at zero temperature can be obtained from the impedance function as 

ReS (u) 



— = 1 — hm ^ 

rrih i^^o CO 



(4.174) 



so that from the scahng relation (I4.172p for the memory function we also have a scaling 
relation for the polaron mass 84j ]: 

(a) "^3D («n«) 



3D 

Since the mobility can be obtained from the memory function as [8 

1 ruh , ImS (u) 
— = lim 



(4.175) 



/i 



e <^-s-o oj 



(4.176) 



fulfilment of the scaling relation fl4.172p implies also a scaling relation for the mobility 84| : 



(J'nD [a) = /i3D [anO!) 



(4.177) 



In the important particular case of 2D, the above scaling relations take the form 82l-l84j : 

(4.178) 



^ , , 2^ /Svr 
E2D [a) = -EsB I —a 



Z2D (a; uj) = ZsY) (^^a; w ) , 
(a) (f «) 



{rrib) 



nD 



[nib, 



3D 



3n 

/i2D laj = yU3D \ —a 



(4.179) 
(4.180) 
(4.181) 
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2. Check of the scaling relation for the path integral Monte Carlo result for the polaron free 
energy 



The fulfilment of the PD-scaling relation 8J] is checked for the path integral Monte Carlo 
results for the polaron free energy. 

The path integral Monte Carlo results of Ref. 43j for the polaron free energy in 3D and in 



2D are given for a few values of temperature and for some selected values of a. For a check 
of the scaling relation, the values of the polaron free energy at /3 = 10 are taken from Ref. 



43| in 3D (Table I, for 4 values of a) and in 2D (Table II, for 2 values of a) and plotted in 



Fig- [m upper panel, with squares and open circles, correspondingly. 



The PD-scaling relation for the polaron ground-state energy as derived in Ref. 8^ reads: 



E2D [a) = -E3D 



37ra; 



(4.182) 



In Fig. [m lower panel, the available data for the free energy from Ref [43| are plotted in 
the following form inspired by the l.h.s. and the r.h.s parts of Eq. (1): F2D {<y) (squares) 
and (^^)(open triangles). As follows from the figure, the path integral Monte Carlo 

results for the polaron free energy in 2D and 3D very closely follow the PD-scaling relation 



of the form given by Eq. ( 4-182 ): 



F2D (a) = -F^D 



Sira 



(4.183) 
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FIG. 11: Upper panel: The values of the polaron free energy in 3D (squares) and 2D (open 
circles) obtained by Ciuchi'2001 [43] for f3 = 10. The data for F^u (a) are interpolated using a 
polynomial fit to the available four points (dotted line) . Lower panel: Demonstration of the PD- 
scaling cf. PDT987: the values of the polaron free energy in 2D obtained by Ciuchi'2001 [43| for 



/3 = 10 (squares) are very close to the PD-scaled according to PD'1987 [84| values of the polaron 
free energy in 3D from Ciuchi'2001 for /3 = 10 (open triangles). The data for ^F^d {^^) 
interpolated using a polynomial fit to the available four points (solid line). 



Appendix 1. Weak coupling: LLP approach 

Inspired by the work of Tomonaga on quantum electrodynamics (Q. E. D.), Lee, Low and 
Pines (L.L.P.) j29| derived fll.22p and m* = mfe(l + a/6) from a canonical transformation 
formulation, which establishes f ll.22p as a variational upper bound for the ground-state 
energy. 
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The wave equation corresponding to the Frohhch Hamiltonian f ll.l6p is 



We shall take advantage of the fact that the total momentum of the system 



(4.184) 



Pop = ^ hkala^ + p 



(4.185) 



(where p = —ikV is the momentum of the electron) is a constant of motion because it 
commutes with the Hamiltonian fll.l6p 
Indeed, 

[p, H] = [p, 5^(\4ake^'^-^ + r;4e-'^")] = J^^V.a^ [P' ^^'''l + [P> ^-^^^'l ) 

k k 

= 5^;ik(rfcake^''-^-V;4e-*^-^); 



J2 hkala^, J^i^k'ak'e"^'-" + V^^le-'^'-' 

k k' 



-ikr\ 



J2hk\Vk 



flk 



4ak, 4 



-5^/lk(l4ake^'^'--l^; 



[Pop,i^] 



0. 



(4.186) 



Because of the commutation (14.1861) . the operators H and Pop have a common set of basis 



functions: = and Pop$ = P^- 



It is possible to transform to a representation in which Pop becomes a "c" number P, 
and in which the Hamiltonian no longer contains the electron coordinates. The unitary 
(canonical) transformation required is $ = 5*1^, where 



5*1 = exp 



^(P - ^ /ikaj^ak)-r 



(4.187) 
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Derivation of the transformations of the operators. 





'p-5i = 


= exp 


. k 


= exp 


k 


= exp 


k 


= P - 


hkalttk + p, 



p exp 



^(P - ^ hkalak)-r 



{—ihS/) exp 



(P-^;ika|,ak)T 



(P - Ek ^kaj^ak) exp [i(P - ^ka|,ak)-r 
+ exp [|(P - Ek^kai«k)-r] (-iW) 



(4.188) 



op 



exp 



= exp 



-|(P - ^fikaj,ak)-r 



^kajj-Ck + p j exp 



/ikajjOk exp 



/ikaj^ak + P - ^kaj^Ck + p = P + p, 



-(P-X^kaj,ak)T 

k 

(4.189) 
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= exp 



ttk exp 



P — /lka|^ak)-r 

k 

— ?ikaj^ak-r 

k 

exp |^ika|^ak-rj Okexp |^— ikaj^ak-r 



exp 



flk exp 



(P — /ikaj^ak)-r 
^J^nkatak-r 



ikaj^ttk-rj ctk^^ ^— ikaj^ak-r j 

ikaj^Ok-rj o"*^*^ ^— ikaj^ak-rj 
ik ■ raj^Okj ^— (— ik ■ r)" Ok (^dkOk 
ik • raj^akj ^ _^—[—ik ■ r ^flkflk + lj]"ak 
ik ■ raj^Okj exp ^—ik ■ r(a|^ak + 1) 



see(*) 



= exp 
= exp 
= exp 
= exp 
= exp 
= Ok exp (— ik ■ r) . 

Here the property was used: 

flk (aia-kj = (aiak + 1 j Ok- 
It is evident for n = O.For n = 1 it is demonstrated as follows: 

QkOk = ^k^k + 1 =^ Ok^kak = (aj^ctk + l)ak; 

then for n > 2 the validity of Q is straightforwardly demonstrated by induction. 
Finally, 



(4.190) 



(*) 



al_ — > Si ^al.Si = [Si ^a^Si]'' = aj^exp (ik ■ r) . 
Using Kim . (1471891) . (ICTOD and flCTT]) . one arrives at 



(4.191) 



H — y H = Si ^HSi 



(p - Xlk ^kaj^ttk 



2mh 



k 



hujLoalak + Y,iVkak + V^al), (4.192) 
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where p is set 0. ^ The wave equation fl4.184p takes the form 



HSiij = ESiij ^HiIj = EiIj. (4.193) 

Our aim is to calculate for a given momentum P the lowest eigenvalue E{P) of the 
Hamiltonian f l4.192p . For the low-lying energy levels of the electron E{P) may be well 
represented by the first two terms of a power series expansion in : E[P) = Eq + P"^ /2mp + 
O(P^), where mp is the effective mass of the polaron. 

The canonical transformation fl4.187p formally eliminates the electron operators from the 
Hamiltonian. LLP use further a variational method of calculation. The trial wave function 
is chosen as 



(4.194) 



where ipo is the eigenstate of the unperturbed Hamiltonian with no phonons present (vacuum 
state). Specifically, ipo is defined by 



OkV'o = 0, (^o,V^o) 
and the second canonical transformation: 



S2 = exp 



(4.195) 



(4.196) 



where /k are treated as variational functions and will be chosen to minimize the energy. 
The physical significance of Eq. (14.196^ is that it "dresses" the electron with the virtual 
phonon field, which describes the polarization. Viewed as a unitary transformation, 5*2 is a 
displacement operator on and aj^ : 



^ Transformation of the equation Pop$ = P$ leads to Si^PopSitp = P-^-At the same time, applying Eq. 
(|4.189p . we obtain Si^PopSi = P + p. Setting the gauge pip = eliminates the operator p from the 
problem. 



88 



exp 



- ak/k) 



ttk exp 



exp [-(4/k - flk/k) 

«k + C^k, (4/k ~ Ctk/k 
flk + /k, 



ttk exp 
1 

' 2 



5^(4/k - flk/k' 

. k 

see(**) 



[flk/k - flk/k) 



Ok, (4/k - Ctk/k) , (4/k - «k/k 



+ ... 



(4.197) 



5-Mk^2 = 4 + /^ 



Here the relation was used 



exp [~V] a exp [V]=a+ [a, V] + ^ [[a, V],V] + ^^ [[[a, V],V],V] + 
Further we seek to minimize the expression for the energy, 



(4.198) 



In virtue of ( 14.1971) and f l4. 19811 . we obtain: 



(4.199) 



Sn ^HS2 



P-Ek^k a^+Z* (ak + /k 



2mh 



+ J2 (ai + f:) (ak + h) + Yl («k + /k) + v: (4 + 



J* - Ek Kkalak) - Ek |/k| - Ek ^k 4/k + Ok/, 



+ ^ ;ja;LO (4«k + |/k|^ + 4/k + flk/k 

k 

+ $^k(«k + /k) + v^fc*(4 + /i 



k 

Ho + iJi, 
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where 



Ho 



|2l2 



hjjj 



k 



■ P 



5^ [Vuk + ^feVk 



2m/, 



nx^LO \- 7; ^ > k /k 



nk-P /i^p ;^2]^ >r-,;,. 
/iwlo \- 7; ^ / k /k 



2mb mb 



(4.200) 



k,k' 



h^k ■ k^ 

2mh 



jakflk'/k/k' + 24ak'/k/k' + 44'/k/k' 



E/i^k • k' I I t t /. 

— <j aj^akflk'/k' + ak'CtkCtkJk' 



k,k' 



2mb 



Using ( I4.195p . we obtain from (I4.199|) that 

p2+ rv. /iki/kl'l' V- 
E = iJo = 2^, + E [^^^ + ^^*/k] 

+ X] l/kl^ S ^Wlo 

k 

We minimize fl4.20ip by imposing 



^k ■ P fi^e 
nib 



5E_ 



6E 



5/. 



This results in 



Vk + f^{ hwi^o + - — + — 

TTib Znib nib 



E^k'l/k 



■k 



(4.201) 



(4.202) 



and the appropriate complex conjugate equation for /k. Upon comparing f l4.202p and f l4.2 



we see that the linear terms in and Ok are identically zero if (14.2020 is satisfied, and 
hence that Hq is diagonal in a representation in which aj^Ok is diagonal. So, the variational 
calculation by LLP is equivalent to the use of (I4.200p as the total Hamiltonian provided 
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satisfies fl4.202p . An estimate of tlie accuracy of the LLP variational procedure was obtained 
by an estimate of the effect of ifi using a perturbation theory. 

Now we evaluate the energy of the ground state of the system, which is given by Eq. 
(I4.2U1I) with satisfying Eq. fl4.202p . The only preferred direction in this problem is P. 
Therefore one may introduce the parameter rj defined as 



(4.203) 



Then Eq. K202^ leads to 



hwi^o — (1 -1]) + 



rrib 



2mh 



(4.204) 



and we obtain the following implicit equation for rj: 



nujLo — (1 -11) + 



2mh 



^ fd'khk(^y^-^(^-] 

{27rf J \ k ) V \2mbUJi^o) 



hjjj 



hk-P 



h^k^ 



LO 



nib ' " 

Let us introduce spherical coordinates with a polar axis along P and denote x = cos(k P): 



riP 



LO 



(- 



h 



27T^ \2mbOJho 



j 2tt / dxx I 



dkk 



hkPx, , h^k^ 
Jtwlo (1 -V) + 



nib 



2mb 



ah f h \ 2 



. . 2-K I dxx / dkk 

2n^ \2nibU)i,oJ J-i Jo 



' ^ hkPx ,^ , ^2^2 

\-2- ^ (1 -r/) + 



2mbhwi,o 



2nibhwi,o 



Further, we introduce the parameter 

Q = 

and a new variable 



{2nibt\uj-Loy''^ 

hk 



(4.205) 



{2nibfiwi,o) 



1/2- 
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This gives 

ah f h 
V = — 

TT \2mbUJLO 



J dxx J dn^K + qx) j '[k^ + {1 — q 



TT P 

\l/2 /.I 



2 2^1 2 



a (2mb^Lo) 



2 2^l2 



g (2mb^Lo) 

TT P 



1/2 .1 L 



-1 I +qx 



2[k2 + (1-£,2^.2)] 



^ arctan 



2(l-g2^2)[^2 + (l_q2^.2)] I 2(l-g2^2)3/2 \^j-^_^2^.2]l/2 



a {2mi,nuJi^of^'^ (1 qxn q^x^ qx • / \ i 

dxx < - + — „ + TTT- + -J- „ arcsm (gx) ^ . 



TT P J_i [2 4(1 - g2a;2)3/2 2(l-g2a;2) 2(1 - g2a;2)3/2 

' dx 



P 4 i_i(l -g2a;2)3/2 



TT 

,"1 2 



4^^-^V-l (l-g^x2)3/2 ^" 



a — yl ~ Q'^ arcsin (g) 

2 g^Y^l — g2 



(2m,,;kJLo)^''^ 



(1-ry))' V/ 



1 - 



— arcsin (g) 



9 — , — arcsin (gj 



So, we have arrived at the equation 



•,(l-lf = ^ (^^^) I ^= - arc=in (,) | , (4.206) 



or equivalently, using the definition fl4.205p . 



« / ^ arcsin (g) 1 . (4.207) 



1-7] 2g3 \ ^/l^ 
Using Eqs. ( I4.203P and f l4.204p . we can simphfy the energy (I4.20ip as follows: 



E = ^^^^^-2y 

r^LO-^(l-r/) + ef 

+ / . ^ ^ ^LO h 



■^(/^LO-^(l-r^) + S^)' 
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{1 + V')P' 2^ \V,f 



1 2 / -D k2u 



+ y]7 — 72 r^LO {'^-ri + ii) + ^ 



(i + /r)p2 



2m, z_.;^,o-lf(l-^) + lS 



k ^LO-ff(l-r7) + S 



^(^LO-^(l-ry) + |g) 



2m, V^^V r^LO-^(l-ry) + e^ 
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The sum over k in Eq. fl4.208p is calculated as follows: 



LO 



V 



V 



dk 



\ V ) 



h 



2mt,aJLO 



IN 2 

4 



fc2 2 
3AI Wlo 



47ra 
1^ 



dk 



dk- 



^LO-^(l-^) + S 



1 



For the calculation of this integral, we can use the auxiliary identity 

• 1 



1 
ab 



1 



t2 • 



Setting 



[ax + b{l — x)y 



(4.209) 



h 



h 



b = e 



we find 



jf2 



7r2 



rribUJloa 



l-ry) + 



2m, i, 



h 

2mbUJi^o 
h \ 



dx / (ik- 



2mbUJ'Lo 

h 



1 



(ix / dk- 



[x (F - 2^^(1 - r/) + + (1 - x) P]' 

1 



dx I dk 



(p_2^f (l-r^)x+ , 

1 



-X 



((k-f (1 - r7)x)' + 2™x - f (1 - Tifx^)' 



dx I dk- 



(F + 2™^x - g(l - r/)2x2)^ 



As long as P^/ (2mf,) is sufficiently small so that no spontaneous emission can occur (roughly 
PV (2"^6) < /ic^Lo), the quantity 

2mbULo 



A = 
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is supposed to be positive for < a; < 1. Therefore, we can use the integral 



what gives 



/^LO-^(l-r^) + 



2 Jo . L (i-'?)'^V 2 



We change the variable x = t^, what gives 



^ 1 , ^ 1 .2 



J X — (f-x"^ Jo Wl — qH'^ Q 



dx = 2 — dt = — arcsing, 



and hence 

> T7-^ ToTT = arcsmg. 4.210) 

As a result, the energy (14.2081) is expressed in a closed form 

E= (1-^) ^arcsing. (4.211) 

Further, we expand the r.h.s. of Eq. fl4.207p to the second order in powers of P (or, what 
is the same, in powers of g) using the relation 

^ arcsin (g) = Ig^ + O (g^) (4.212) 



what results in 



T] a 
1 — rj 2g^ 



+ O (?=) 



Solving this equation for rj, we find 

l + a/6 ylmbfkjOi^o J 



(4.213) 
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We also apply the expansion in powers of g up to ~ to the energy fl4.21ip : 
(l - rf) - ahuLo - -ahuLoq^ + It^loO (g"^) 



2mft ^ 6 

p2 

= -a/k^LO + Y^;^ (1 - r/) ((6 + a) r/ - a + 6) + HuloO (g^) 

p2 

= -ahwLo + 77- 77T + ^LoO (g^) . 

2mb{l + a/6) 

Finally, we have arrived at the LLP polaron energy 
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Appendix 2. Expansion in Stieltjes continuous fractions 



5^ 



In this derivation it is shown that the approximation used in the evaluation of the function, 



which determines the polaron mass [see Eqs. (40) and (Bl) from Ref. 52| ] 

g{k,z) = / dre'^^exp [-^0(0)] exp [k^Cir)] (4.215) 

J —oo 



with 



C{r) = \j:-,\h'f'^ (4.216) 

is equivalent to an expansion in a continued fraction limited to the first step. Moreover, it 
is proved that the choice of the coefficients of the continued fraction can be justified by a 
variational principle, at least when z is real and positive. 

Expanding the last exponential of Eq. ( I4.215p in a power series leads to 



g{k,z) 



.0 oo . X n 

^ y. fcffc|...t;iA,i'i/fei'...iA„i' 



X exp [i(7fci + 7fc2 + •••TfcJ^] 

= -I exp [-PC(O)] / (4.217) 



n=0 

^ klkl..kl\f,f\hf...\f,f 1 

X 



"'I ) ■ ■ ■ 

The multiple sum over the fc's is in fact an integral with 3n variables. It is possible to change 
the variables in that one of the new variables is 

Xn = 7fci + 7fc2 + ••• + 7fc„- (4.218) 

Then the multiple sum which appears in the last term of Eq. fl4.217p is of the following 
type: 

oo 

J{z) = [ ^^dxn, (4.219) 

J Xn + Z 
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where 



L{Xn) = 

is the resuh of the integration over the n — 1 other variables. An expansion of integrals of 
the type (14.2191) into Stieltjes continued fractions is known to give good results when z is 
real and not located on the cut of J{z), i.e., when 

z > -nu. (4.220) 
The first nontrivial step in the continued fraction expansion is 



J(z) = (4.221) 
^ ' ai + z ^ ^ 



with 



oo 



ao = / L{xn)dxn, (4.222) 



nu] 



oo 



a,^"-^ . (4.223) 



A variational principle can be established, which gives a rather strong argument in favour 
of the approximation (14.2211) . Let us introduce a variational parameter x writing 

oo 

TILU 

Performing two steps of the division, this relation becomes 



oo 

J{z) = J L{Xn)dXn 

oo 

- T—^ I - x)L{xn)dXn + K{z, x) (4.225) 
{z + xY J 

nuj 
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with 



oo 

K{z, x) = / ^^^^^^^rfx.. (4.226) 



The approximation consists of neglecting the term K{z,x) in Eq. fl4.225p . As this term 
is positive [cf. fl4.220p ]. the best approximation is obtained when it is minimum. Therefore 
let us use the freedom in the choice of x to minimize the expression f l4.226p . 



dK{z,x) ^K{z,x) 



dx z + X 

CO 

-2-4-Ti / ^^^^^^^4^^^" = 0' (4-227) 

[Z + X)'^J Xn + Z 



which gives 



or 



oo 

2^^ / + ^^L{x^)dxn = (4.228) 

{Z + Xj-'J \ Z + X J Xn + Z 



oo 

1 



[Z + X) 



- / (x„ - x)L{xn)dxn = 0. (4.229) 



This provides us with the best value of the variational parameter 



XyiLi^Xfi^dXyi 

x=-^ , (4.230) 



which is ai [cf. Eq. ( KTH^ ' . 

With this value of x and neglecting K{z,x), the expression fl4.225p of the calculated 
quantity J{z) becomes 

oo 

J{z) = [ L{xn)dxn = J{z) = (4.231) 

z + X J ai + z 



which is the first step (I4.22ip of a Stieltjes continued fraction. 
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To prove that this value of x gives a minimum of K{z,x), let us calculate the second 
derivative 



d'^K{z,x) 6 



oo 



dx"^ (z + x) 



^Xn x^ Zji^x^^dxri 



nu) 
oo 



+ T-^-T^ I L{xn)dxn. (4.232) 



^3 
nuj 



{z + x) 

Now the parameter x is replaced by its expression fl4.230p . The relation fl4.232p becomes 

d'^K{z,x) 2 



oo 



j L{xn)dxn, (4.233) 



\3 

nu) 



dx^ {z + ai] 

which is positive of z ^ —nu, since it follows from relation f l4.230p that ai > nu. 

Our approximation is related to that used by Feynman which is based on the following 
inequality: 



e""^) ^ e-'<^\ (4.234) 
where the brackets denote the expectation value of the random variable x. For instance, 



L{x)e '^^dx 

{e-'^) = , (4.235) 

L{x)ds 



where L(x) is the non-normalized probability density of x. The Laplace transform of Eq. 
(0341) gives 

nOO POO 

/ e~'' {e-'"") ds ^ / e-^"e-^<">(is, (4.236) 
Jo Jo 

which after integration becomes 



oo 



r/ \ I L(x)ds 

^^"^ dx ^ = (4.237) 



la X + z ^ (x) + z ai + z 
The last inequality shows the relation with our procedure. 
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Part II 

Many polarons 



V. OPTICAL CONDUCTIVITY OF AN INTERACTING MANY-POLARON GAS 
A. Kubo formula for the optical conductivity of the many-polaron gas 



The derivations in the present section are based on Ref. 51|- The Hamiltonian of a 
system of N interacting continuum polarons is given by: 

N 2 

j=i q 

N 2 ^ ^ 1 

q i=i °° i=i £(^j)=i ' ^' 

where rj, pj represent the position and momentum of the constituent electrons (or holes) 
with band mass m^; &q,&q denote the creation and annihilation operators for longitudinal 
optical (LO) phonons with wave vector q and frequency wlo! ^ describes the amplitude 
of the interaction between the electrons and the phonons; and e is the elementary electron 
charge. This Hamiltonian can be subdivided into the following parts: 

H = He + i^e-e + Hph + He-ph (5.2) 

where 

N 2 

^-E^ (5.3) 

i=i 

is the kinetic energy of electrons, 

2 N N 



E E ^ (5.4) 



is the potential energy of the Coulomb electron-electron interaction, 

^Vh = ^LO&q&q (5.5) 
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is the Hamlitonian of phonons, and 



N 



(5.6) 
q i=i 

is the Hamiltonian of the electron-phonon interaction. Further on, we use the second quan- 
tization formahsm for electrons, in which the terms He, He-e and He-ph are 



k,o- 



2mb 



He-e = ^ E^^l E 4+q,a4'-q,a'«k',<x'ak,a = ^ E ' ^^^-^^ 
q^^O k,<T q^O 

k' ,cr' 

q 

where : ... : is the symbol of the normal product of operators, 

47re^ 
''''^^^ 

is the Fourier component of the Coulomb potential, and 



TV 

Pq = E 



'^k+q,o-'^k,a- 



(5.7) 
(5.8) 

(5.9) 
(5.10) 
(5.11) 



k,(T 



is the Fourier component of the electron density. 

The ground state energy of the many-polaron Hamiltonian (15. ip has been studied by L. 
Lemmens, J. T. Devreese and F. Brosens (LDB) 89], for weak and intermediate strength of 
the electron-phonon coupling. They introduce a variational wave function: 



IV'ldb) = U 



(0) 

el I ' 



(5.12) 



V'el'^ ) represents the ground-state many-body wave function for the electron (or hole) 



where 

system and |0) is the phonon vacuum, ?7 is a many-body unitary operator which determines 
a canonical transformation for a fermion gas interacting with a boson field: 



I j=l q 



(5.13) 



In Ref. [89|, this canonical transformation was used to establish a many- fermion theory. 
The /q were determined variationally 89|] resulting in 



/q 



2^2 ' 



hrq 



(5.14) 



LO 



2m65(q) 
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for a system with total momentum P = Y2j Pj = 0. In this expression, 5'(q) represents the 
static structure factor of the constituent interacting many electron or hole system : 

N N 

^^(q) = ( ^^^e^"'^"'^)- (5-15) 




The angular brackets (...) represent the expectation value with respect to the ground state. 

The many-polaron optical conductivity is the response of the current-density, in the 
system described by the Hamiltonian (15. ip . to an applied electric field (along the x-axis) with 
frequency u. This applied electric field introduces a perturbation term in the Hamiltonian 
(15 . 1 p . which couples the vector potential of the incident electromagnetic field to the current- 
density. Within linear response theory, the optical conductivity can be expressed through 
the Kubo formula as a current- current correlation function: 

""^"^^^'v^^Vh^ J, ^'^'^[Ut).JM])dt. (5.16) 

In this expression, V is the volume of the system, and is the x-component of the current 
operator J, which is related to the momentum operators of the charge carriers: 

with q the charge of the charge carriers (+e for holes, — e for electrons) and P the total 
momentum operator of the charge carriers. The real part of the optical conductivity at 
temperature zero, which is proportional to the optical absorption coefficient, can be written 
as a function of the total momentum operator of the charge carriers as follows : 

Re a{uj) = -^^Re{ r e^-* ([P,(t), P,(0)]) dt] . (5.18) 



103 



B. Force-force correlation function 



Let us integrate over time in (15.181) twice by parts as follows: 

) 



1 



ico — S ^ 



U=0 



dt 



d 



iuj — 5 Jq 
1 



iuj — 5 
1 

iu — 5 Jq 
1 



ioj — 6 Jq 

1 

iu — 6 



dt 

) 

dt 
dt 
dt 



dt 

Awt—5t 



Acjt—St 



h 



Aiot—St 



Auit—St 



F,(0),e--^P,e 



t=0 



dt 



d 



Awt—St 



iujt—St 



1 ^ 


iu — 








iuj — 












i 


/ 



P.(0),e-f^(^[/7,P.])e^^ 



([P,,P,])+ / dt 

Jo \ L 

([P., PJ) + / rft ( P, (0) , e-f ^P,. (0) e 



Au}t—St 



h 



([P,., PJ) + / dt( [e^^P,. (0) e-^^, P, (0) 



Aujt—St 



iu}t—5t 



Awt—St 



{lo + id) I Jo 
where F = I [P, P] is the operator of the force applied to the center of mass of the electrons. 

Since both P^. and P^ are hermitian operators, the average ([Pr, Px]) is imaginary. Hence, 
for 7^ 0, this term does not give a contribution into Re a {u) . As a result, integrating by 
parts twice, the real part of the optical conductivity of the many-polaron system is written 
with a force-force correlation function: 

1 



Recrfw) 



Re|^"e^"*([P,(t),P,.(0)])dt 



(5.19) 



The force operator is determined as 



P. = ^ [H, PJ = ^[He + Pe-e + Pp. + Pe-p., Px] ■ 
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Further, we use the commutators: 



k' 

^ ^ j «k+q,a«k,aak/,^ak',a + a 
k' \ -^k'.a^k+Q.a'^k'.aai 



4+q,a'^k',<7«k,aak' 

-^k',a^k+q,f7'2k',<7ak, 



k' 
k' 



^/,^Ctk',<7ak+a.a"k,a 



Wq,o-'^k,c 
'k',k+qak',o-'^k,<7) 



^+q,t7«k,a - 4',k+qak+q,<^Ctk,,7) 



= a 



tk+q,^ak,a XI ("^kk' - 4',k+q) = -^?xak+q,^ak,<x, 

k' 



Hence, [i/e, P.] = 0, [H^-e, Px] = 0, 

[ffe-ph, Px] = Y (^A [Pq> + [P-q, Pj) 

q 

= X]^^ (Ki&qPq - Kj*^qP-q) , 

q 

So, the commutator of the Hamihonian (15. ip with the total momentum operator of the 
charge carriers leads to the expression for the force 

F = E q i^^^^P^ - KKp-^) ■ (5-20) 
q 

This result for the force operator clarifies the significance of using the force-force correlation 
function rather than the momentum-momentum correlation function. The operator product 
Fx{t)Fx{0) is proportional to |l^kP; the charge carrier - phonon interaction strength. This 
will be a distinct advantage for any expansion of the final result in the charge carrier - phonon 
interaction strength, since one power of |Vkp is factored out beforehand. Substituting (I5.20p 
into (I5.19p . the real part of the optical conductivity then takes the form: 

Vhw^ml J, ^^^^^ 

X {[[V^b^ (t) + ^4^;^ (t)] pq (t) , (y_q,6_q, + V*,b+) p_q,] ) . (5.21) 

Up to this point, no approximations other than linear response theory have been made. 
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C. Canonical transformation 



The expectation value appearing in the right hand side of expression fl5.2ip for the real 
part of the optical conductivity is calculated now with respect to the LDB many-polaron 
wave function (I5.12p . The unitary operator ( I5.13P can be written as 

f/=eXpJ]AqPq, = /q&q - /^q&lq, (5.22) 

q 

The transformed Hamiltonian (I5.2p is denoted as 

H = U-^HU. (5.23) 

The momentum operator of an electron p^, the operator of the total momentum of elec- 
trons P and the phonon creation and annihilation operators are transformed by the unitary 
transformation (15.220 as follows: 

U-^PjU = + J] hciA^e'^^ ^'i, (5.24) 
q 

U-^VU = P + 5^ ^q^qPq, (5.25) 
q 

U-XU = 6q - /*p_q, U-XU = b+- /qPq. (5.26) 



As a result, after the transformation fl5.22p . the Hamiltonian takes the form (see Ref. 89|): 

H = He + He-e + Hph + He-ph + Hm + -f/ppe, (5.27) 

where the terms are 

^e-e = ^ • P^P-^ •' ^q = ^q + 2 [^Lo |/q|' - Vq/* - V*f^) , (5.28) 

q^O 



He-ph = [(Vq - hWLoU) KP^ + K - ^Lo/q) ^qP-q] 

q 

+ ^ E E + 2k • q) 4+q,.«k,., (5.29) 

q kjCr 

Hn = nY, (^lo |/q|' - Vj; - V;U) , (n^ E<<^«k,J , (5.30) 

q V k,a / 
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(5.31) 



qq' 



The exact expression for the real part of the conductivity (15.211) after the replacement of 
|\E'o) by I^^ldb) = U \ (f)) is transformed to the approximate one 



Re a (u) 



1 e' 



Re 



q,q' 



u 



;(0) 



1 e' 



Re / dte'^^'-^'^q^q'^ 



q,q' 



; (0) 



Re / dte'^'-'' 



Vhuj'^ ml 



q.q' 



{V^ (feq - /qP-q) + K^q (^^q " /-qP-q) ) Pqe"^^, 



- /!q,pq,) + V*, (6J - A,pq,)) p_q,] I 0)1 V^^)) 



So, we have arrived at the expression 

Re a (u) 



1 p2 /-oo _ 

Re / cite'"*-^* 



mt 



e*"" (K, (&q - /qP-q) + K!, (fel, " /-qP-q)) Pqe'^ ^, 
(V-q' (fe-q' - Tq'Pq') + V*, (foj - f^, p^,)) p_q,] | 0> | 



,(0) 



Since p^P-^ = P-^Pq, and Kj/q = ^-q/*q, the terms proportional to P-^p^ vanish after the 
summation over q: 



'i^^^f*oiP-^P^ = - X] ^^^q/qP-qPq = 0- 



(5.32) 



107 



Hence we obtain the real part of the optical conductivity in the form 



Re a ioj) 



Re 



dte 



ILUt—St 



Yl ^^^^ 



q,q' 



(y_q,6_q,+\/*,6+)p_q,]l0>l^: 



,(0) 



(5.33) 



Introducing the factor 



I (0) 



(5.34) 



the optical conductivity can be written as 



Re a iu) 



(5.35) 



In the case of a weak electron-phonon coupling, we can neglect in the exponent e~^^ of 
f l5.33p the terms He-ph and Hppe [i- e., the renormalized Hamiltonian of the electron-phonon 
interaction f l5.28p and f l5.3ip ]. Namely, we replace H in Eq. fl5.33p by the Hamiltonian 



In this case, we find 



Ho = He + He-e + Hph + H^- 



J(q, q') = (^1?) U e^^o (K,6, + 1^461 J p^e"^^". 



(5.36) 



y-q'&-q'+K;'&j)p-q']k>l^. 



/ (0) 



2i|lqp5qq' Im 



;(0) 



e^^»*/Vq&qe-'^"*/V-qfe: 



The time- dependent phonon operators are 



/(O) 
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so that we have 



:r(q,qO = 2z|Vq|2V Im 



2i\Vq^\ 5qq' Im 



;(0) 



(0) 







where He = He + He-e + Hm. 

Taking the expectation value with respect to the phonon vacuum, we find 



J(q,q') = 2z|Vqp5qq, Im 



I (0) 



(5.37) 



The optical conductivity (I5.33P then takes the form: 

2e2 



Re a {uj) 



Im 



dte 



iujt—St 



X Im 



el 



;(0) 



(5.38) 



For an isotropic electron-phonon system, in 3D can be replaced by | (g^ + + g^) = |g 
what gives us the result 

e 



ReasD (w) 



2 -2 



POO 

Vg^lVql^Im / dte''^^-^Hm[e-'^^''^F{q,t 



(5.39) 



where the two-point correlation function is 



it fjr i^ fj 

eft^^PqC -""^P-q 



The same derivation for the 2D case, provides the expression 
where A is the surface of the 2D system. 



(5.40) 



(5.41) 



D. Dynamic structure factor 



To find the formula for the real part of the optical conductivity in its final form, we 
introduce the standard expression for the dynamic structure factor of the system of charge 
carriers interacting through a Coulomb potential. 



iq.(rj(t)-rKO)) 



(5.42) 
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The dynamic structure factor is expressed in terms of the two-point correlation function as 
follows: 



1 

2N 



(0) 



-oo 
oo 



— /" F (q, t) e^^'dt = —F (q, oj) 



S{q,u) = ^Fiq,u) 



where F (q, u) is the Fourier image of F (q, t): 



F{q,u)= j F (q, t) e'^'dt. 

— oo 

The function F (q, t) obeys the following property: 



/ (0) 
/ (0) 



it fj it fj 



it TT 

Pqe" V-q 



where Eq is the eigenvalue of the Hamiltonian Hg 



(0) 



I (0) 



Herefrom, we find that 



F*(q,t) = e-*^o/^(o) 



He 



PqC'* "P-q 



/ (0) 



It rr Z-T 



^(q.-t). 



From f l5.45p . for the function 

B (q, t) = Im [e-'^'Loti? (q^ t) 
the following equality is derived: 

B (q, -t) = Im [e'^^otp 

= Im [e'^^ot;p* (q^ 

= - Im [e-'^'^o'F (q, t)] = -5 (q, t) 



(5.43) 



(5.44) 



(5.45) 



(5.46) 



110 



B{q, -t) = -B (q,t). 

The integral in Eq. fl5.39p 

noo POO 

Im / dte''^^-^^ Im [e-'^^o^F (q, t)] = Im / dte'^^-^^B (q, t) 
Jo Jo 



(5.47) 



is then transformed as follows: 

Im / dte"^^-^^B (q, t) = — 
Jo 2z 

1 

~ 2i 
1 



dte''^i-^'B{q,t) - I dte-''^'-''B{q,t) 



Uo 



iujt—5t 



dte''^'^^'B{q,t) - I dte''^'+''B{q,-t) 



Uo 



iujt+St 



-oo 




dte'^'-^^'B (q, t) 



rfte*^*"^*fi(q,t) - 

— j rfte*'^*-''l*l5(q,t) 

2^ J -oo 2i 

— / rfte*'^*-''!*! [e-^'^L°*F (q, t) - e''^'^°'F (q, -t)] 

4 J —oo 



We can show that, as far as 



is the ground state, the integral j'^^dte^^* ^l*lF(q, — t) 



for positive a; is equal to zero. Let | V'lr^^j the total basis set of the eigenfunctions of 
the Hamiltonian He. Using these functions we expand F(q, t): 



'(0) 



it fj a fj 



(0) 



rfte^--^l*lF(q,-t)=E|(^lr)|p_,|^, 

n 

= El('''er'ip-<,i'/', 



i?r, — En 



2n6 \uj + 



0, 



because for > 0, w + is never equal to zero. 

So, rewriting expression f l5.39p with the dynamic structure factor of the electron (or hole) 
gas results in: 



oo 
oo 



1 /"°° 
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and we obtain 



Recr3D(w) 

Re(T2D(w) 



no 



ml 



^g^|Vq|^S'(q,u; - Wlo), 



no \ - 



2hu^ ml ^ 



g^|Fq|^S'(q,w - Wlo) 



(5.48a) 
(5.48b) 



where 



no 



N/V in 3D, 
N/A in 2D 
is the density of charge carriers. 

For an isotropic medium, the dynamic structure factor does not depend on the direction 
of q, so that S'(q, w) = S{q,u), where q = |q|. Let us simphfy the expression f l5.48ap using 
exphcitly the amphtudes of the Frohhch electron-phonon interaction. The modulus squared 
of the Frohhch electron-phonon interaction amplitude is given by 



IK 



(^Iwlo)^ 47ra; f h 



V \2mbUJi^o 



in 3D 
in 2D, 



q A \2mbULo 
where a is the (dimensionless) Frohlich coupling constant determining t 
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(5.49) 



le coupling strength 



between the charge carriers and the longitudinal optical phonons 
respectively, the sums over q is transformed to the integrals as follows: 

V 



8^. In 3D and 2D, 



q 

2D:$:... 



{2ny 

A 
(2^ 



dq. . . 
(iq . . . 



4 



Rea3D(w) 



no V 
3hw^ ml (27r) 



3 / dq^q^ 



iq \ V 



1/2 



h 



2mbUjLo 



1/4 



S{q,uj- Wlo) 



Re(T3D(w) 



noe^ 2a ^jjIq f h 



ml 37r u'^ \2mhUi^o 



1/2 



q^S{q,uj- Ui,o)dq. 



(5.51) 



In the same way, we transform Re (120(1^): 

, , noe^afkolc, f ^ f°° 2r,^ n, 

RecT2D(w) = — 2-77 / q S{q,u -u-Lo)dq. 

m^ 2 u-^ \2mbULoJ Jo 



112 



Using the Feynman units {h = 1, nib = 1, i^lo = 1), Reo"(c<;) is 

^V^a 1 



Rea3D(w) = noe^— q'^S{q,uj - l)dq, 

Reo-2D(w) = ^oe^— 7=-? / q^S{q,UJ - l)dq. 
zy Jo 



2^2 

From these expressions, it is clear that the scahng relation 

Recr2D(w,a) = Re cr3D(c<;, — a) 



(5.52) 
(5.53) 



(5.54) 



which holds for the one-polaron case introduced in ref. 



83 



84| . is also valid for the many- 



polaron case if the corresponding 2D or 3D dynamic structure factor is used. 



1. Calculation of the dynamic structure factor using the retarded Green's functions 



The dynamic structure factor S (q, u) is expressed through the two-point correlation 
function by Eq. f l5.43p . The correlation function F(q, w) can be found using the retarded 
Green's function of the density operators 



G^(q,t) = -.e(t) 



(5.55) 



where G (t) is the step function. Let us consider the more general case of a finite temperature, 



(q,t) = -iQ{t) 



where the average is 



(...) 



/3 



Tr {e-f^^o) ' ksT' 
The Fourier image (q, u) of the retarded Green's function G^ (q, t) is 



(5.56) 



(5.57) 



/oo 
G^ (q, t) e'^'dt 
-oo 



P-qe^^Vqe-^'^'' ) ) e'-'dt 
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The imaginary part of then is 

J 

J 



^ -'-OO ^ II 7-0O ^ 



In the second integral here, we replace t by (i' + 

Tr (|e-^^V-qe^^Vqe"^^°) 



/ 

J — c 



^n. / _V -^e-*'dt'. 



As far as the integral over t converges (i. e., ^e^^^Pqe ft^°p_qy tends to zero at \t\ — >■ 
we can shift the integration contour to the real axis, what gives us the result 



Tr (e-l^"°p-qe'i"°pqe-'i"° 



e'^^'dt 



.00 Trfet^°-^^Vqe-^'^V-q') 
J-00 Tr (e-/3^o) 

t 



L 



— e 



(ef^Vqe-^^'V-q) e'^'dt. 

■00 

Herefrom, we find that 

ImG^ (q,a;) = (l - e"'^^"^) H (e^^Vqe-t^o^_ \ ^iu^t^^ 

Im G« (q, a;) = (1 - e"^'^) F (q, a;) . 
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So, the equation follows from the analytical properties of the Green's functions: 

2ImG^ (q,w) 



F{q,u) 



(5.58) 



The formula (15.581) is related to arbitrary temperatures. In the zero-temperature limit 
(/3 —7- oo), the factor (l — e~^^'^) (I5.58P turns into the Heavicide step function 6 [u), what 
leads to the formula 



^(q,w)lr= 



-26 (w)Im (q,w) 
4 



lT=0 



(5.59) 



S{q,u)\^^, = --e{co)lm G«(q,u;) 



\t=o (5-60) 

The retarded Green function is related to the dielectric function of the electron gas by the 
following equation: 



G^(q,a;) = - 



1 



Within the random phase approximation (RPA), following 



(5.61) 



70|, the expression for (q, u) 



IS 



{q,u) = [1 - v^P (q, u)]-' hP (q, u) , 
where the polarization function P (q, w) is (see, e. g., p. 434 of [tO]) 

1 \ ^ /k+q,cr ~ /k,(T 

h 2^ 



P{q,uj) 



5 -> +0 



(5.62) 



(5.63) 



with the Fermi distribution function of electrons /k,o-. 

For a finite temperature, the explicit analytic expression for the imaginary part of the 
structure factor P3D (q, w) is obtained (see |70|). 



Im P3D (q, uj) = 



Vml j^l + exp{/3[/i-E(+)(g,a;)]} 



27rn4/3g 1 + exp {/3 [/i - P(-) (g, w)] } 

(^^sr 



4^ 



(5.64) 



with the chemical potential fi. The real part of the structure factor is obtained using the 
Kramers-Kronig dispersion relation: 



Re P (q, u) 



V 



1 

u' — u 



ImP (q, u') doj' . 



(5.65) 



115 



Analytical expressions for both real and imaginary parts of P (q, u) can be written down 
for the zero temperature (see [70|), 



Re P3D {q, uj) 



Vnih 



+ 



TFq^ \ 2m, 







hq^ 
2m^ 


hkpq 




In 






hq'^ 


, hkpq 









ImPsD {q,uj) 



Vrrih 



h2q2 + 2mt 



+2kFq 

2mi J 





In 


1 hq^ 
2 '"6 


hkpq 




1 hq'^ 


hkpq 









2m(, 



u-2 



2mi, J 



(5.66) 



1 /3 

where kF = {37r'^N/V) is the Fermi wave number. 

After substituting into Eq. f l5.58p the retarded Green's function f l5.62p in terms of the 
polarization function we arrive at the formula 

2h ^ P{q,uj) 



F{q,uj) 



Im 



S (q,w) 
'S'(q,w)lr=o 



h 



ImP (q,a;) 



AT (1 - e-/^'--) [1 -v^ReP (q, u)f + K Im P (q, iu)f ' 

h ImP(q,a;) 

N ' [l-t;qReP(q,a;)]' + KlmP(q,a;)]' 



(5.67) 
(5.68) 



With this dynamic structure factor, the optical conductivity (15.521) (in the Feynman units) 
takes the form 



Re 0-30(0;) = -e' 



V2a 1 
3nVu^ 



X 



ImPao {q,uj - 1) 



'0 [1 - t;gReP3D (g,w - 1)] + [wglmPjD (g,w - 1)] 
Correspondingly, in the 2D case we obtain the expression 

a 1 



jq'^dq. 



(5.69) 



2V2V 



:e (uj - 1) 



X 



ImPaD (g,w - 1) 



[1 - t;gReP2D (g,w - 1)]^ + [t;gImP2D (g,w - 1)]^ 



q^dq. 



(5.70) 
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2. Plasmon-phonon contribution 



The RPA dynamic structure factor for the electron (or hole) system can be separated in 
two parts, one related to continuum excitations of the electrons (or holes) 5*00111; and one 
related to the undamped plasmon branch: 



SRPA{q,uj) = Api{q)S {u - Wpi (g)) + Scont{q,uj) 



(5.71) 



where Upi (g) is the wave number dependent plasmon frequency and Api is the strength of 
the undamped plasmon branch. 

In Eqs. fl5.69p . fl5.70p . the contribution of the continuum excitations corresponds to the 
region {q,uj) where lmP{q,uj) 7^ 0. The contribution related to the undamped plasmons is 
provided by a region of (g,^;) , where the equations 



Im P{q,uj) = 
1-VgReP (g, u) = 

are fulfilled simultaneously. Using fl5.6ip . we find that Eqs. f l5.73p are equivalent to 



(5.72) 



Im 



e{q,uj) 



0, Re 



e{q,uj) 



0. 



(5.73) 



In the region where Im P (q,a;) = 0, the expression rj- 



lmP{q,uj) 



-Vg Re P{q,uj)y + [vg Im P(g,w)] 



^ IS propor- 



tional to the delta function, which gives a finite contribution to the memory function after 
the integration over g: 



ImP (g,^) 



[l-VgReP{q,oof +[vglmP{q,uj)f 



5{l-VqReP{q,uj)). 

TTVn 



lmP{q,ui)=0 



(5.74) 



Using Eq. f l5.74p . the coefficients Ap\{q) in Eq. fl5.7ip can be expressed in terms of the 
polarization function P (g, u) as follows: 

ImP {q,uj) 



[l-VqReP{q,ujf+[vglmP{q,uj)f 



^^.NIjReP(g, 



00 



ImP(g,w)=0 
6{UJ- Wpl (g)) 



'^='^pl{<?) 
4 
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^pl (?) 



(5.75) 



='^pl(9) 



The contribution derived from the undamped plasmon branch A^\{q)5 {u — ujp\{q)) is 
denoted in Ref. jsil as the 'plasmon-phonon' contribution. The physical process related to 
this contribution is the emission of both a phonon and a plasmon in the scattering process. 



E. Comparison to the infrared spectrum of Nd2_a;Cej;Cu02-j/ 

Calvani and collaborators have performed doping-dependent measurements of the infrared 
absorption spectra of the high-T^ material Nd2_xCe^Cu02-j; (NCCO). The region of the 
spectrum examined by these authors (50-10000 cm~^) is very rich in absorption features: 
they observe is a "Drude-like" component at the lowest frequencies, and a set of sharp 
absorption peaks related to phonons and infrared active modes (up to about 1000 cm~^) 
possibly associated to small (Holstein) polarons. Three distinct absorption bands can be 
distinguished: the 'd-band' (around 1000 cm^^), the Mid-Infrared band (MIR, around 5000 
cm~^) and the Charge- Transfer band (around 10^ cm"^). Of all these features, the d- 
band and, at a higher temperatures, the Drude-like component have (hypothetically) been 



90|. 



associated with large polaron optical absorption 

For the lowest levels of Ce doping, the d-band can be most clearly distinguished from 
the other features. The experimental optical absorption spectrum (up to 3000 cm~^) of 
Nd2Cu02-5 (5 < 0.004), obtained by Calvani and co-workers JooJ, is shown in Fig. [T2] 
(shaded area) together with the theoretical curve obtained by the present method (full, 
bold curve) and, for reference, the one-polaron optical absorption result (dotted curve). At 
lower frequencies (600-1000 cm~^) a marked difference between the single polaron optical 
absorption and the many-polaron result is manifest. The experimental d-band can be clearly 
identified, rising in intensity at about 600 cm~^, peaking around 1000 cm^^, and then 
decreasing in intensity above that frequency. At a density of no = 1.5 x 10^'' cm~^, we found 
a remarkable agreement between our theoretical predictions and the experimental curve. 
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FIG. 12: The infrared absorption of Nd2Cu02-<5 < 0.004) is shown as a function of frequency, 



90] is represented by the 



up to 3000 cm^^. The experimental results of Calvani and co-workers 
thin black curve and by the shaded area. The so-called 'd-band' rises in intensity around 600 cm~^ 
and increases in intensity up to a maximum around 1000 cm^^. The dotted curve shows the single 
polaron result. The full black curve represents the theoretical results obtained in the present work 
for the interacting many-polaron gas with no = 1.5 x 10^^ cm~^, a = 2.1 and = 0.5 mg. (From 
Ref. [51].) 



F. Experimental data on the optical absorption in manganites: interpretation in 
terms of a many-polaron response 



In Refs. 



91 



92| . the experimental results on the optical spectroscopy of La2/3Sri/3Mn03 
(LSMO) and La2/3Cai/3Mn03 (LCMO) thin films in the mid-infrared frequency region are 



presented. The optical conductivity spectra of LCMO films are interpreted in |91l . |92| | in 
terms of the optical response of small polarons, while the optical conductivity spectra of 
LSMO films are explained using the large-polaron picture (see Fig. [T^ . 



The real part of the optical conductivity Re a (u) is expressed in Ref. |92| by the formula 

2 {huof 



Re (j{uj) 



arir. 



'S^''\q,u-uJo)dq, 



(5.76a) 



where a is the electro n-phonon coupling constant, Up is the polaron density, m is the electron 
band mass, Uq is the LO-phonon frequency, and S^^\q,u) is the dynamic structure factor, 
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FIG. 13: Comparison of the low-temperature MIR optical conductivity to a (u) from various model 



calculations: the solid line refers to the weak-coupling approach of Tempere and Devreese [51 1 
modified for an on-site Hubbard interaction, the dashed line is the result of the phenomenological 
approach for self-trapped large polarons by Emin |93| , the dotted curve is the weak-coupling single- 
polaron result [gO^. (From Ref. [oi].) 



determined in Ref. |92| through the dielectric function of an electron gas e (q, cu): 



CO 



In Ref. 



Up 47re^ 



Im 



(5.77) 



51| , the other definition for the dynamic structure factor is used, which is equivalent 



to that given by fl5.77p (see Ref. 70|) with the factor (the number of electrons) 

+ 00 

5(™)(q,a;) = i /" (^ei|Pq(t)p-q(0)|¥.ei)e^^*rft = iV^(^)(q,a;). 



(5.78) 



Here, (pei denotes the ground state of the electron subsystem (without the electron-phonon 
interaction), j9q = X^jLi e^^'^^ is the Fourier component of the electron density. 

In Ref. j5l|, the dynamic structure factor is calculated within the random-phase ap- 
proximation (RPA) taking into account the Coulomb interaction between electrons with the 
Fourier component of the Coulomb potential 



(5.79) 



where ^oo is the high-frequency dielectric constant of the crystal. In Refs. [9l|,l9^, the dy- 



namic structure factor is calculated taking into account the local Hubbard electron-electron 
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interaction instead of the Coulomb interaction. The local Hubbard interaction is used in 
the small-polaron formalism and describes the potential energy of two electrons on one and 
the same site (see, e. g., Refs. |9J, |95|). In its simplest form the Hubbard interaction is (see 
Eq. (1) of Ref. [Mi) 

VH = Uj2^,^n,^, (5.80) 

i 

where U is the coupling constant of the Hubbard interaction, rii^ is the electron occupation 
number for the i-th site. 



In Refs. 



91 



92l |. there are no details of the calculation using the interaction term (15.801) . 



The following procedure can be supposed. The transition from the summation over the 
lattice sites to the integral over the crystal volume V is performed taking into account the 
normalization condition 

where the density n^(i) (r) is to be determined through As far as the lattice cell volume 

VLq ^ V , the integral (r) dr can be written as the sum over the lattice sites: 

/ «ta) (r) rfr = fio V ^ta) (^i) , 
Jv , 

where {rj} are the vectors of the lattice. Therefore, from the equality 



we find that 

The potential f l5.80p is then transformed from the sum over sites to the integral: 

Vff = Tii^nnQo = I Uno6{r- r ) (r) (r) drdr'. 

"0 "~ Jv 



(5.81) 



Consequently, in the continuum approach the Hubbard model is described by the 6-like 
interparticle potential UQq6 (r — r'). 

This development of the approach 5l| performed in Refs. 9l|, seems to be con- 
tradictory by the following reason. For a many-polaron system, both the electron-phonon 
and electron-electron interactions are provided by the electrostatic potentials. Therefore, 
it would be consistent to consider them both within one and the same approach. Namely, 
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the Coulomb electron-electron interaction with the potential (15.791) is relevant for large and 
small polarons with the Frohlich electron-phonon interaction, while the Hubbard electron- 
electron interaction is relevant for small Holstein polarons. Nevertheless, as recognized in 
Ref. jojl, this model "reproduces the observed shape of t he p olaron peak quite convincingly" 
and provides a better agreement with the experiment 
approach [93| and the one-polaron theory 



91 



92l | than the phenomenological 
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VI. INTERACTING POLARONS IN A QUANTUM DOT 



A. The partition function and the free energy of a many-polaron system 



We consider a system of electrons with mutual Coulomb repulsion and interacting 
with the lattice vibrations following Ref. 96|]. The system is assumed to be confined by 



a parabolic potential characterized by the frequency parameter Qq. The total number of 
electrons is represented as = Yla ^o"' where No- is the number of electrons with the 
spin projection cr = ±1/2. The electron 3D (2D) coordinates are denoted by x^ o- with 
j = 1, ■ ■ ■ , Nfj. The bulk phonons (characterized by 3D wave vectors k and frequencies u^) 
are described by the complex coordinates Qi^, which possess the property [41 1 



Qk — Q-k- 



(6.1) 



The full set of the electron and phonon coordinates are denoted by x = {xj^o-} and Q = {Qk} ■ 
Throughout the present treatment, the Euclidean time variable r = it is used, where t is 
the real time variable. In this representation the Lagrangian of the system is 



L (k, Q- X, Qj = Le (^, x) - Vc (x) - Ub (x) + ^ [Q, + L^-ph (x, Q) , (6.2) 
where L^, (x, x) is the Lagrangian of an electron with band mass mf, in a quantum dot: 

L4i,x) = - >; > ; :^ (A|„ + nijxy , a^:^, (6.3) 



VLq is the confinement frequency, Vb (x) is the potential of a background charge (supposed to 
be static and uniformly distributed with the charge density en^ in a sphere of a radius i?), 

N 



VAr) 



cr j = l 

350 



QR2_ 2 p3- 

Q{r<R) — + Q(R<r) — 

2 r 



(6.5) 



where is the static dielectric constant of a crystal, Vc (x) is the potential energy of 
the electron-electron Coulomb repulsion in the medium with the high-frequency dielectric 

constant ^oo: 

a,a'=±l/2 j = l 1=1 °° ' ^''^ ' 
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Lph (Q, Q*', Q, Q*^ is the Lagrangian of free phonons: 

(O, q) = -\ 5^(QkQk + c^^QLQk), (6.7) 

k 

Further, Le-ph (x, Q,Q*) is the Lagrangian of the electron-phonon interaction: 

1/2 



Le-ph (x, Q) = - ^ f j VkQ-kPk, 



(6.8) 



k 

where pk is the Fourier transform of the electron density operator: 

P^= H E^^'"'^'^- (6-9) 

cr=±l/2 j = l 

Vk is the amplitude of the electron-phonon interaction. Here, we consider electrons inter- 
acting with the long-wavelength longitudinal optical (LO) phonons with a dispersionless 



frequency uj^ = wlo, for which the amplitude W is [24 1 



Vi^^fH^V'V^)'". (6,10) 

Q \ V J \mbUjLoJ 

where a is the electron-phonon coupling constant and V is the volume of the crystal. 

We consider a canonical ensemble, where the numbers N„ are fixed. The partition func- 
tion Z {{Ncr} , P) of the system can be expressed as a path integral over the electron and 
phonon coordinates: 

Z({N,}.P)=Z [,l;'X,,.'. I * /r^"'"* / '''5^%«We-*"*'l. (6.11) 

where S [x (r) , Q (r)] is the "action" functional: 

S [x (r) , Q (r)] =-^J^ L Q- x, (Qj dr. (6.12) 

The parameter /3 = 1/ {ksT) is inversely proportional to the temperature T. In order to 
take the Fermi-Dirac statistics into account, the integral over the electron paths {x(r)} 
in Eq. fl6.1ip contains a sum over all permutations P of the electrons with the same spin 
projection, and denotes the parity of a permutation P. 

The action functional (16.121) is quadratic in the phonon coordinates Q. Therefore, the path 
integral over the phonon variables in Z {{N„} , (3) can he calculated analytically following 
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Ref. j4l|. Let us describe this path integration in detail. First, we introduce the real 
phonon coordinates through the real and imaginary parts of the complex phonon coordinates 
Q'^^ = ReQk, Qk = Ii^Qk- According to the symmetry property (16.11) . they obey the 
equalities 

Q-k = Q'k> Q\ = -Ql (6.13) 
gk = < ^ (6.14) 

[ V2Qk, K < 0. 

In this representation, the sum over phonon coordinates IQ^f is transformed in the 



following way using the symmetry property (16.131) : 

E iQki' = E [(Q'^f + (^k 



k)^ 



k k 

(fc:.>0) (fe:,<0) 

= E E ^k = E^k- 

k k k 

(fcx>0) (fc^<0) 

Therefore, the phonon Lagrangian (16. 7p with the real phonon coordinates is 

^p/^ = -^E(^k + ^kgk)- (6.15) 

k 

The Lagrangian of the electron-phonon interaction (16. 8p with the real phonon coordinates 
is transformed in the following way using (I6.13P : 

1/2 



'^e-,H = - E ^kPk (g'k - ^Qk) 

= - E (^) ^^^p^ + ^-^p-^^ Q'^ 



k 

(fc:.>0) 



E (^)'^'(^kPk-vikP-k)Q;:. 



k 

(fc:.<0) 

Let us introduce the real forces: 

TA'-^)y^p^^y-^p-^h^^^^^ (6.16) 

173 (¥) (^kPk - V-kP-k) , fc.>0. 
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This gives us the Lagrangian of the electron-phonon interaction in terms of the real forces 
and real phonon coordinates: 



■'e—ph 



7kgk- 



(6.17) 



So, the sum of the Lagrangians of phonons and of the electron-phonon interaction is ex- 
pressed through ordinary real oscillator variables: 



Lph + L 



e—ph 



45Z(^k + ^k?k + 7kgk)- 



(6.18) 



The path integration for each phonon mode with the coordinate gk is performed indepen- 



dently as described in Sec. 2 of Ref. 



Dgk (t) exp 



4l| and gives the result 

1 [^^ 1 
~ io 2 ^ ^^^^ ^ ^^^^^ 



2sinh(^) 



'^^^^^coshic^kdr-r'l- ^2)] 



cjk sinh {PhuJk/2) 



7k [r) 7k [r 



where the exponential is the influence functional of a driven oscillator { 411], Eq- (3.43)}. 
Therefore, the path integral over all phonon modes is 



//•{9k} 
d{q^} / 

n 



^{?k(^)}exp 



1 /■'"'^ 1 

^ X] 2 ^ ^^^^ ^ ^^^''^ 



2 sinh (^) 



X exp <i ^ / dr 



hp 



cosh [wkdr-r'l -hl3/2)] 
Wk sinh (/3/k^k/2) 



7k j 7k 



Here, the product Yl^ • • • is the partition function of free phonons, and the exponential is 
the influence functional of the phonon subsystem on the electron subsystem. This influence 
functional results from the above described elimination of the phonon coordinates and is 
usually written down as e~*, where $ is 



hi3 

■Tn' * 



hl3 



k 



cosh [tJk ( I 



r — T 



cjk sinh {l3huj]i/2) 



-7k [r) 7k [T 



(6.19) 
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The sum over the phonon wave vectors k can be simphfied as follows: 
^cosh[uk{\T - t'\- hl3/2)] 



k 



- y 



sinh {f3huj]^/2) 
cosh[wk(k-r'| -hl3/2)] 



h ^ sinh (/3/iCi;k/2) 

(fc:.>0) 



X [FkPk (r) + V^-kP-k (r)] [VkPk (r') + V^-kP-k (r')] 
1 ^ cosh[wk(k-r'| -n/3/2)] 

~ n ^ sinh (/3/kJk/2) 

(fc^<o) 

X [l^kPk (r) - r-kP-k (r)] [^kPk (r') - V^-kP-k )] 
^ 2 ^ cosh[a;k(|r-rV;i/3/2)] 
h ^ sinh (/3;iWk/2) ^ 

X [Pk (r) p-k (r') + p-k (r) Pk (r')] 

2.^ coshK(|r-rVn/3/2)] , 
= 12-^ sinh(/3;^k/2) l^*^' Pic(r)p-k(r). 

Herefrom, we find that 

V 2^ ■/o io sinh(^) ^''^ ki ; i j 

As a result, the partition function of the electron-phonon system (16.111) factorizes into a 
product 

z({iv.},/^) = z,({iv.},/^)n ,,,..h('r^,/2) 

of the partition function of free phonons with a partition function Zp {{Ncr} , /3) of interacting 
polarons, which is a path integral over the electron coordinates only: 

ZA{N,},P)^^j^^^^ld^fD^(r)e~^''-^^^\ (6.22) 

The functional 

Sp [x (r)] = -1 ^'"^ [L, (r) , X (r)) - Vc (x (r))] + $ [x (r)] (6.23) 

describes the phonon-induced retarded interaction between the electrons, including the re- 
tarded self-interaction of each electron. 
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Using (16. 3p and (16. 6p we write down Sp [x (r)] explicitly: 



5,[x(r)] 



h 



o-.ct' j = l «=1 

(i,<x)7^{/,<x') 

hp/2)] 



2^2 



cr J=l 

(ir / (ir 



2£, 



,cosh [wlo (I 



T — T 











sinli {/3hLUi,o/2) 



(6.24) 



The free energy of a system of interacting polarons Fp {{N^} is related to their par- 
tition function (I6.22p by the equation: 



({iVj, /?) = -- In Z,({iVj,/3). 



(6.25) 



At present no method is known to calculate the non-gaussian path integral (16.2211 analyt- 
ically. For distinguishable particles, the Jensen- Feynman variational principle |4l| provides 
a convenient approximation technique. It yields a lower bound to the partition function, 
and hence an upper bound to the free energy. 

It can be shown [97] that the path-integral approach to the many-body problem for a 
fixed number of identical particles can be formulated as a Feynman-Kac functional on a 
state space for indistinguishable particles, by imposing an ordering on the configuration 
space and by the introduction of a set of boundary conditions at the boundaries of this state 
space. The resulting variational inequality for identical particles takes the same form as the 
Jensen- Feynman variational principle: 



F < F 



-^0 + {Sp - So)g^ , 



(6.26) 
(6.27) 



where Sq is a model action with the corresponding free energy Fq. The angular brackets 
mean a weighted average over the paths 

Ep / Z^x (.) (.) e-"[^Wl 

((•))5o = ^ ' . . .p. _ . . . .. . ■ (6-28) 



fd5tf!'''D5t{T) e-^o[x(r): 



In the zero-temperature limit, the polaron ground-state energy 



= }^ Fp 

p—>-oo 



(6.29) 
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obeys the inequality following from f l6.26p with f l6.27p : 
with 

.0 - n 



E^^r = El^ + \im ^- {Sp - So)s,j , (6.30) 
E° = lim Fq. (6.31) 



B. Model system 

We consider a model system consisting of electrons with coordinates x = {xj^o-} and 
Nf fictitious particles with coordinates y = {y^} in a harmonic confinement potential with 



elastic interparticle interactions as studied in Refs. 96|, |98| . The Lagrangian of this model 
system takes the form 



,2 



Lm (x, y; X, y) = ^ ^ {±1^ + r^'x^ J + " ^'.-')' 

cr j = l o-.cr' j = l 1=1 

- ? E (y] + "/y]) - 1 E E E - yf ■ ^^■^^) 

j = l a j = l 1=1 

The frequencies fi, w, the mass of a fictitious particle m^, and the force constant k 
are variational parameters. Clearly, this Lagrangian is symmetric with respect to electron 
permutations. Performing the path integral over the coordinates of the fictitious particles 
in the same way as described above for phonons, the partition function Zq ({A^o-} , /?) of the 
model system of interacting polarons becomes a path integral over the electron coordinates: 

ZMN,),ff) = Y,Jj^J^J <ix£''m(r)e-*I>«l (6.33) 

with the action functional Sq [x (r)] given by 

s» (^)i = i /"' E E T [*- + * 

a j = l 

-his N^, 2 



1 f rribUj 
Jo „ - ,_i 



cr,a' j = l 1=1 

hi3 hp 

f 



/• , /• , ,cosh[fi.(|r-r'| -^1/3/2)1 , , , 

Arrifhrif J J smh (/3/ifi//2) ^ ^ ^ ^ ^ ^ 
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where X is the center-of-mass coordinate of the electrons, 

X=iv5ZE^^- (6.35) 

(7 j = l 

1. Analytical calculation of the model partition function 

The partition function Zq {{N„} , /3) [Eq. f l6.33p ] for the model system of interacting 
polarons can be expressed in terms of the partition function Zm ({^o-} , iV/, /?) of the model 
system of interacting electrons and fictitious particles with the Lagrangian Lm [Eq- (I6.32p 
as follows: 

where Zf {Nf,Wf,(3) is the partition function of fictitious particles, 

1 

(2 sinh iPhwf) 

with the frequency 



Zf (^/' = ....... .DNr (6.37) 



Wf = ^Qj + kN/rrif (6.38) 

and D=3(2) for 3D(2D) systems. The partition function Zm {{N^^} , Nf, (3) is the path 
integral for both the electrons and the fictitious particles: 

//.px f fy 

dx J D5t{T) J dy J Dy(r)e-^«[^(")'y(")l (6.39) 



with the "action" functional 

-hp 



Sm [x (r) , y (r)] = ~^ J (^^ y) dr, 



(6.40) 



where the Lagrangian is given by Eq. fl6.32p . 

Let us consider an auxiliary "ghost" subsystem with the Lagrangian 

(X„ Y„ X„ Y,) = [±l + ^^Xj) - ^ (y^^ + w}Y^) (6.41) 

with two frequencies w and Wf, where w is given by 



w = - Nu^ + kNf/rrib. (6.42) 
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The partition function Zg of this subsystem 

= 1 d^g J dYg I DXg (r) I DYg (r) exp {-5, [X, (r) , Y, (r)]} , (6.43) 



with the "action" functional 



S»[X,M.Y,(.)l^-i/L,(x,.X,„Y„Y,).. (6.44) 







is calculated in the standard way, because its Lagrangian fl6.4ip has a simple oscillator form. 
Consequently, the partition function Zg is 



Zg = ^ n- (6-45) 

[2sinh(^)]"f2sinh^^^^" 



The product ZqZm of the two partition functions Zg and Zm ({^o-} , P) is a path 
integral in the state space of electrons, Nf fictitious particles and two "ghost" particles 
with the coordinate vectors and Y^. The Lagrangian Lm of this system is a sum of Lm 
and Lg, 

Lm (x, y, Xg, Y^; x, y, X^, Y^^ = Lm (x, y; x, y) + Lg (^X^, Y^, X^, Y^^ . (6.46) 

The "ghost" subsystem is introduced because the center-of-mass coordinates in Lm can be 
explicitly separated much more transparently than in Lm- This separation is realized by the 
linear transformation of coordinates, 

= X;. , + X - Xg, ^^^^^^ 

y,. = y;., + Y-Y„ 

where X and Y are the center-of-mass coordinate vectors of the electrons and of the fictitious 
particles, correspondingly: 

cr j = l J j = l 

Before the transformation f l6.47p . the independent variables are (x, y,X^,Yg) , with the 
center-of-mass coordinates X and Y determined by Eq. f l6.48p . When applying the trans- 
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formation fl6.47p to the centers of mass fl6.48p . we find that 



X = ^ E E (-;-.^ + X - X,) = ^ x;,^ + x - x„ (6.49) 

(6.50) 



a j = l a j = l 



Nj Nf 

^=W, S + Y - Y,) = -i- ^ y;, + Y - Y, 

As seen from Eqs. fl6.49p . fl6.50p . after the transformation fl6.47p the independent variables 
are (x',y',X, Y) , while the coordinates (Xg, Y^) obey the equations 

x.=^EE4.. Y,=-i-Ey;- (6-51) 

o- j = l J j = l 

In order to find the explicit form of the Lagrangian (16.460 after the transformation (16.471) . 
we use the following relations for the quadratic sums of coordinates: 

E E 4. = E E (4 J ' + (X= - XJ) , 5] yj = 5; (y^ ) ^ + AT, (Y^ - Y,=) , 

cr j = l a j = l j = l j = l 

E E E - = E E (4.) " - 

a,a' j=l 1=1 a j = l 

Nf Nf Nf 



r2 

■ 9' 

i=l 1=1 3=1 



Na N„ Nf 

E E E (-i^. - yf = ^/ E E + 'vE W)' + 

a j=l 1=1 a j = l j = l 

A^A^/(x2 + Y2-2X-Y-XJ-Y^). (6.52) 

The substitution of Eq. (I6.48P into Eq. (I6.46P then results in the following 3 terms: 

Lm (x , p', X, Y; x', y', X, y) = (k' , x') + L^^. (p' , y') + Lc (x, X; Y, y) , (6.53) 

where ^x',x'j and L^f ^y',y'j are Lagrangians of non-interacting identical oscillators 
with the frequencies w and Wf, respectively, 



X') = E E + K ' (6-54) 

<T=±l/2 j = l 
Nf 
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The Lagrangian Lc ^X, X; Y, describes the combined motion of the centers-of-mass of 
the electrons and of the fictitious particles, 



with 



Lc (X,X;Y,Y 



5 ^5 -^5 



Y' + wjY' + fcATATyX ■ Y, (6.56) 



n = + kNf/mb. (6.57) 

The Lagrangian fl6.56p is reduced to a diagonal quadratic form in the coordinates and 
the velocities by a unitary transformation for two interacting oscillators using the following 
replacement of variables: 

1 



with the coefficients 



ai 



X 



X 
Y 

l + X 



1 



(air + a2R) , 



/ 



1/2 



02 



1-x 



1/2 



1 1/2 ' 



(Ci^-nfj +472 

The eigenfrequencies of the center-of-mass subsystem are then given by the expression 



NN 



f 



mi,mf 



(6.58) 

(6.59) 
(6.60) 



\ 



\ 



n'^-n^f] +47 



+ 472 



(6.61) 



As a result, four independent frequencies Qi, Q2, w and Wf appear in the problem. Three of 
them (fii, ^^2, w) are the eigenfrequencies of the model system. Vti is the frequency of the 
relative motion of the center of mass of the electrons with respect to the center of mass of 
the fictitious particles; VL2 is the frequency related to the center of mass of the model system 
as a whole; w is the frequency of the relative motion of the electrons with respect to their 
center of mass. The parameter is an analog of the second variational parameter w of the 
one-polaron Feynman model. Further, the Lagrangian (16.561) takes the form 



(6.62) 
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leading to the partition function corresponding to the combined motion of the centers-of- 
mass of the electrons and of the fictitious particles 

1 1 



1 D ■ 



(6.63) 



[2sinh [2sinh(»)]^ 

Taking into account Eqs. fl6.45p and fl6.63p . we obtain finally the partition function of 
the model system for interacting polarons 

D 



sinh (^) sinh 
sinh (») sinh (») 



Zp{{N^},w,(3) 



Here 



Zf ({ATJ ,w,f3) = Zp {Ni/2, w, /3) Zf w, /?) 



(6.64) 



(6.65) 



is the partition function of = Ni/2 + N_i/2 non-interacting fermions in a parabolic confine- 
ment potential with the frequency w. The analytical expressions for the partition function 
of spin-polarized fermions Zp {N^j, w, (3) were derived in Ref. |99| . 



C. Variational functional 

In order to obtain an upper bound to the free energy E^ar, we substitute the model action 
functional f l6.34p into the right-hand side of the variational inequality f l6.26p and consider 
the limit /3 — )■ cxd: 



E.ar ({iVj) 



N 



Ef ({ivj ,w) + ^{nl-n' + Nu') ( J2 ^] (0) 



So 



(X^ (0)),^ + {U, (x)),„ + ^ [g (q, 0| {iVj , /3 ^ oo) - AT] 



hp np 

k^N^Nf f f ,cosh[^],(|r-r'| -^/3/2)] s ^ , 
+ hm ^ / dr / dr' ^ , , . (X (r) • X (r')) 





ri/3 hp 



sinh {/3hnf/2) 



So 



r /"w f , cosh[a;Lo(|r-/|-fl/3/2)] ^^ 'UMin^ 
/ I sinh(/3;^,o/2) S{ci,r-r\{N^},P). 







(6.66) 
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Here, 'Kf{N,w) is the energy of non-interacting fermions in a parabolic confinement 
potential with the confinement frequency w, 

(7=±l/2 I n=0 ^ ^ 

+ {N^-N0(^L^ + l^Y (6.67) 

where a is the spin of an electron, Lg. is the lower partly filled or empty level for iVo- electrons 
with the spin projection a. The first term in the curly brackets of Eq. fl6.94p (the upper 
line) is the number of electrons at fully filled energy levels, while the second term (square 
brackets) is the number of electrons at the next upper level (which can be empty or filled 
partially). The energy levels of a 3D oscillator are degenerate, so that 



g(n) = -{n + 1) {n + 2) 
is the degeneracy of the n-th energy level. The parameter 



(6.68) 



Nl. = -K + 1) + 2) 

6 



(6.69) 



is the number of electrons at all fully filled levels. The summation in Eq. f l6.67p is performed 
explicitly, what gives us the result 



{L, + if (L, + 2) + (iV. - iV^J (l, + ^ 



(6.70) 



In Eq. f l6.66p . ^ (q, t — r'| {N^} , (3) is the two-point correlation function for the electron 
density operators: 

G (q, r| {iV4 , /3) = (pq (r) (0))^^ . (6.71) 
The averages (X (r) • X (t'))^^ are calculated using the generating function method: 



{X,{t)X,{t')),^ 



(X(r)-X(r')) 



m ^ of cosh [a (|r - a| - n/3/2)] 



€=0, 



So 



2mN 



i=l 



Qi sinh {hf3Q^/2) 



(6.72) 



(6.73) 
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Substituting this expression into Eq. f l6.66p and performing integrations over r and a ana- 
lytically, we obtain the result 

hp hi3 

k'N'Nf r r cosh[nf{\T~T'\-hf3/2)] 

J sinh{mf/2) (XW-X(r)),^ 





3fi7 of 



4 tr^/-^^ 



coth(/3fii/2) coth (/3^]//2) 



r2j r2 J 

and in the zero-temperature limit we have 

,,. cosh[l],(|r-r-|-W2)] . ^ 



(6.74) 



WV^ . . cosh in, fig/2)] (.) . X (r')). 







/TV \ 

The average ( ^ is transformed, using the described above operations with the 

/so 

"ghost" subsystem, 

X,- = x; + X - Xg, (6.76) 
and taking into account the first of equations fl6.52p 

N N 

E-? = EK)' + ^(X=-XD- (6-77) 

i=i i=i 

Consequently, averaging the left-hand side of Eq. fl6.77p on the model action functional Sq, 
one obtains 

N \ / N \ IN 



fxA +iv((x^)^^-(xj)^j. (6.78) 

N \ 

The term ( X] / expressed using the virial theorem through the ground-state energy 
E,p {N, w) of independent 3D fermion oscillators with the frequency w and with the mass 

ij = l / Q 

On; 
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Two other terms in Eq. f l6.78p are [cf. Eq. fl6.73p ]: 



2mbN 



3h j2 co^h (/3^i/2) 



1=1 

3h coth(/3w/2) 



^^^9 2mhN 



w 



So, we obtain 



N 

'^=1 ' So 



2 .2 



Ef{{N^}, w) ^ 3h l^af 1 



(6.80) 



(6.81) 



The averaging of the background-charge potential gives us the result 



1)' + 2 



(T n=0 



k=0 



k\ \n — k / \ 2w 



^fc-l/2 



r] = 6oo/£o, A = 



iFi{k- -; -; ] - iF^ { k - 

^ ^ ' 2' 2' / ^ M 2' 2' 



(6.82) 



i=l 



iV-1 



where /i (ra, (t|/3, A^^-) is the one-particle distribution function of fermions (the distribution 
functions are considered in more details in the next subsection). 

Collecting all terms together, we arrive at the variational functional 



W 



nl + w^ E{w,N) 3 
2^2 h 2 



f) 



+ {Ub {^))s, +Ec + Ee-ph, 

where Ec and E^-ph are the Coulomb and polaron contributions, respectively: 



(6.83) 



E, 



c 



E. 



V2a 



(iq- 



oo 



g{q,0\{K},f3)\ 



/3— 5>oo 



-N 



e—ph 



(6.84) 



j rfq^ j dTexp{-uj^or)g{ci,T\{N^},l3)\p^^. (6.85) 



The correlation function f l6.7ip is calculated analytically in the next subsection. With 
this correlation function, the variational ground-state energy is calculated and minimized 
numerically. 
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D. Two-point correlation functions 



The two-point correlation function f l6.7ip is represented as the following path integral: 



1 ^ f-l)«^ 



Zo({iV4,/3)^iVi/2!iV_i/2! 

xjd5tj_ Dx(r)e-^o[^("Vq(7-)P-q(0). (6.86) 

We observe that Q (q, t| {A^o-} , /3) can be rewritten as an average within the model "action" 
Sm [x (r) , y (r)] of interacting electrons and fictitious particles: 



X 



d5t / Dx(r) I dy I Dy (r)e-'5'>^[^M'y(")] 



xpq(r)p_q(0). (6.87) 

Indeed, one readily derives that the elimination of the fictitious particles in f l6.87p leads 
to fl6.86p . The representation fl6.87p allows one to calculate the correlation function 
^ (q, r| {Ai'o-} , /3) in a much simpler way than through Eq. (16.860 . using the separation 
of the coordinates of the centers of mass of the electrons and of the fictitious particles. This 
separation is performed for the two-point correlation function (I6.87P by the same method as 

it has been done for the partition function ( 16.390 . As a result, one obtains 

(exp[iq- (X (r) -X(a))])q 
g (q, r| {iVj ,P) = g (q, r| {iVj , /3) / ' ' x / m ' ^^.88) 

(exp[«q- (Xg (r) -Xg (a))])^^ 

where g (q, r| {N^j} , /3) is the time-dependent correlation function of N non-interacting elec- 
trons in a parabolic confinement potential with the frequency w, 

~g (q, r| {N,} , (5) = {p^ (r) (0))^^ . (6.89) 

The action functional [x^] is related to the Lagrangian (x, x) [Eq. (I6.54p 

Sw [xr] = ^jL.aj (x, x) dr. (6.90) 



The averages in (I6.88P are calculated using Feynman's method of generating functions 



4l|. Namely, according to ^], the average 



G[/(r)] = (^exp|^y"/(r)x.dr||^ , (6.91) 
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where S^j is the action functional of a one-dimensional harmonic oscillator with the frequency 
u and with the mass m, results in 



G If (r)] = exp {-^fdrf , /°^M^(|r -^1 - /?/2)] ^ ^ ,^^^2, 



Amhhj 



sinh (/3a;/2) 





The diagonalization procedure for the Lagrangian Lc f l6.56p allows us to represent that 
Lagrangian as a sum of Lagrangians of independent harmonic oscillators, what gives the 
following explicit expressions for averages in Eq. fl6.88p : 



(exp[^q-(X(r)-X(a))])5^ 

%2 " 



exp 



Nnih 



2 sinh 

2 



E 

i=l 



^i|T-g-| 

2 



sinh 



a(^./3-|r-(7|) 
2 



l]^sinh(^) 



(exp[zq- (X, (r)-X,(a))])^^ 

%2 sinh sinh 



exp 



w{nP-\T-CT\) 



Nruh 



w 



sinh (^) 



i. T/ie correlation function 5 (q, t| {A'^o-} , /3) 



As seen from the formula (16.891) . g (q, r| {A'o-} , /3) is the time-dependent correlation func- 
tion of non-interacting fermions in a parabolic confinement potential with the frequency 
w. Let us consider first of all a system of identical spin-polarized oscillators with the 
Lagrangian 

N 

^ /',v2 , ,2„2^ 

i=i 

The corresponding Hamiltonian is 

N 



L = -^(x--a;-x 



(6.93) 



/ ^ 9 9 ^ '? 



H = y — H 

i=i ^ 



(6.94) 



N 



,2^2 



(6.95) 



A set of eigenfunctions of the one-particle Hamiltonian h is determined as follows: 



htpn (x) = Eni^n (x) 



(6.96) 
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where 



n 



[ni,n2,n3) 
Sn = Sn = ^ i n + 



Ipn (x) = iprn i^l) V^na (3^2) {x^) 



(6.97) 



n = ni + n2 + ris, 

(fn {x) is the n-th eigenfunction of a one-dimensional oscillator with the frequency u. 

The Hamiltonian f l6.94p can be written down in terms of the annihilation (Sn) and creation 
(a^) operators: 

= ^ SnCi^an = ^nA^n, = a+fln- (6.98) 



The many-particle quantum states in the representation of "occupation numbers" are 
written down as | . . . A^n • • • ) , where A^n is the number of particles in the n-th one-particle 
quantum state. The states | . . . A^n • • • ) are defined as the eigenstates of the operator of the 
number of particles in the n-th state A'^n: 

Nn\...Nn...) = Nn\...Nn...). (6.99) 

Let us determine a set of quantum states with a finite total number of particles 

^Nn = N (6.100) 

n 

as follows: 

1--N^---)\j:^n^=n^\^n,{n^})- (6.101) 

Further on, we use the basis set of quantum states (16.1011) for the derivation of the partition 
function, of the density function and of the two-point correlation function. 

Partition function 

The density matrix of the canonical Hibbs ensemble is 

p = exp(-f3H), (3 



The partition function of this ensemble is the trace of the density matrix on the set of 
quantum states (I6.10ip : 



Z,(/3|Ar) = 5^(M/^,l^„l 



exp -I3H 



N,{Nn} 



exp [ -f3 enN„ 

{Nr.} 



(6.102) 
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This expression can be written down also in the form 



{Nr.} 

where 



Zj{P\N) = J]exp L/3^£„ivJ 6n,j:^n^, (6.103) 



1, J=k 
^0, Jy^k 

is the deha symbol. 

Let us introduce the generating function for the partition function in the same way as in 
Ref. |99i: 

oo oo / \ 

s {13, u) = J2 ^""^i = E E -/5 E^-^" '^^'E.^n 

N=0 N=0 {Nn} V n / 

{Nn} \ n J N=0 



{A^n} \ n / 



S ^) = n ^ E (-/5^n)]^" \ . (6.104) 



Fermions 



For fermions, the number A^n can take only values A^n = and A^n = 1- Hence, for 
fermions (denoted by the index F), we obtain: 



Ep {/3, m) = JJ [1 + M exp (-/?£:„)] • 



Since the n-th level of a 3D oscillator is degenerate with the degeneracy 

(n + 1) (n + 2) 



gin) 



2 



we find that the generating function Ep (/3,m) is given by 

oo 

(/3, u) = l[[l + u exp (-/3£„)]^("^ . (6.105) 



n=0 

Bosons 
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For bosons (denoted by the index B), = 0, 1, ... , oo. The summations over {A^n} in 
Eq. fl6:TO gives: 



2b (/3,w) = n 



n=0 



1 - wexp {-I3en)\ 

The results (I6.104p and fl6.106p prove (for the partition function) the equivalence of the 
path-integral approach for identical particles [99| and of the second-quantization method. 



Integral representation 



Let us use the Fourier representation for the delta symbol: 



2tt 



d9. 



(6.107) 



where ( is an arbitrary constant. Substituting Eq. fl6.107p into Eq. f l6.103p we obtain 

2tt 



Z,(/3|iV)=5^exp(-/3 5^£„iV„|i- /" exp J^iV„-iV| 

{A^n} \n /^o L\n / 



(0 - <) 



d9 



^ f dOexp [-ziV(^-zC)] J^exp ( -/? ^ £„iV„ + ^ ^ iV„ (^ - zC) ] 

^ i^"} V n n / 

^ j de exp {-iN9 - NC) S {(3, e^^+^) ^ 



27r 



Zi{f3\N) = — [ d9exp HnS (/3, e^^+^) - N( - iN9] . (6.108) 
27r J 



The partition function for a finite number of particles can be obtained from the generation 
function also by the inversion formula 100 1 



2m T 2^+1 



1 

2^ 



2n 



[in H(/3,tje''' ) -TV In n] ^-iNO 



(6.109) 
(6.110) 



Let us denote in Eq. f l6.108p : 

C^lnw. (6.111) 

In these notations, Eqs. f l6.108p and f l6.110p are identical. For the numerical calculation, it 
is more convenient to choose in Eq. (I6.107P the interval of the integration over 9 as [— tt, tt] 
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instead of [0, 2n] , what gives: 



Zj{p\N) = — / ^Ni0)de, 

Zn 



(6.112) 



with the function 

^n{0) =exp[\nE{/3,ue'^) -Nlnu-tNO] . (6.113) 

The aforesaid method of derivation of the partition function [Eqs. fl6.107p to fl6.108p 
is heuristically useful, because it allows a simple generalization to spin-mixed systems with 
various polarization distributions. 

The two-point density- density correlation function in the operator formalism is 



^(q,r|{iV4,/3) = (pq (r)/)_q(0)), 
where pq {t) is the density operator in the Heisenberg representation: 

Pq (^) = exp (^^H^ pqexp (-^^) • 
In the "second-quantization" representation, pq (t) is 



n,n' 



(6.114) 



(6.115) 



(6.116) 



After substituting Eq. fl6.116p into fl6.114p . we find that 



h 



n,n' m,m' 

(6.117) 

The operator a^a^'a^a^' has non-zero diagonal matrix elements in the basis of quantum 
states |\l'Ar,{Ar„}) Only in the cases 



n = n I n = m 

or 

m = m' I m = n 



(6.118) 



Hence, the average 



N,{Nn} 



(6.119) 
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is not equal to zero only when the condition fl6.118p is fulfilled. This allows us to write down 
the average (16.1191) as 



1 

(1-6- 



= ^n'n^m'm (1 - ^mn) ^NnNm) + (^m'n^n'm (1 - ^mn) ( ( ^^n ) - ( N^N, 
'^n'n'^m'm'^mn-^n 

Here, the notation is used for the average occupation number A^n^ 



(6.120) 



Z,(/3|iV) = 5^(vl>^,|;v„} 



exp [-I3H)N^ 
exp (-PH 



(6.121) 



^exp I -/3^£„Ar„ j 



(6.122) 



In the same way as Eq. (I6.102p . the average (16.1211) can be written down in the form 



Zi {I3\N) 



Ar„exp {-l3Y^e^,N^i 

{N^,} \ 



m (/3|iV) 



/3Z,(/3|iV) fe„ 



27rZj {I3\N) J exp (/3e„ - C - i^) + 1 



de. 



(6.123) 



Since = ^n, \^n/ depends only on n. 
Using Eq. (I6.11ip . we can write Nn as 

1 } {e) 



Nr. 



-de 



2-kZi{I3\N) J i exp (/3£, - ^^) + 1 

— TT 

1 } exp [In S we^^) - iV In m - iNO] 



— TV 

1 J S T/e^^) exp [-tO (iV - 1) - 

2'kZi{(3\N)J l + Mexp(i^-/3e„) 



(6.124) 
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The averages ( A^n^m / for m 7^ n can be also expressed in terms of the integral repre- 
sentation: 



NnN^ 



Af„,| \ n' / 



Af,{Afn} 



mT^n 



(/3|iV) 
1 



52Z, (/3|iV) 



1 



6' 



2tt 



Zi {/3\N)/3^ 6e^6e 



- I ^ J dOexp [inS e*^+^) - A^C - iNO] 



We obtain the integral representation for the average of the product of operators NnNm for 
m 7^ n: 



(9) 



2ttZi (N) J [exp - C - i^) + 1] [exp -(-16) + !] 

— TV 

Let us introduce the notation 



(19. (6.125) 



1 



(6.126) 



exp {l3e-C-i9) + V 

which formally coincides with the Fermi distribution function of the energy e with the 
"chemical potential" {( + i9) / [3. Using this notation, the averages ( I6.123p and f l6.125p can 
be written down in the form 



Nr. 



1 



Nr. Nr. 



2tiZi{(5\N) 
1 



m^n 2nZi (N) 



f{en,9)<!>N (9) d9, 
f{en,9)f{em,9) $^ (9) d9. 



(6.127) 



(6.128) 



We can develop the aforesaid procedure for the average of a product of any number of 
operators iVn^iVnj . . . iVn^, where all quantum numbers ni,n2,...,n^ are different. The 
result is: 



N^,N^,...N, 



1 



n,^n, 27rZ7 (N) 



f {en„9) f {e^,,9) ...f (£„^, 9) {9) d9. (6.129) 
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It should be emphasized, that all expressions above [including Eq. fl6.129p ] are derived for 
a canonical Hibbs ensemble (i. e., for a fixed number of particles) and for both closed-shell 
and open-shell systems. 

Let us substitute the average (16.1201] into the correlation function g (q, r| {N^} , 



g (q, r| {K} ,^) = J2I1 (^"■')nn' (^""■')„.„.' 



n,n' m,m' 



m'm 



~g (q, r| {iVj , /3) = (e^'^-^)„^ (iV„iV, 



n,m 



The matrix elements fe"''') has the following form 

\ / nm 



(6.130) 



where (m le*"^^] ra) is the matrix element of a one-dimensional oscillator with the frequency 



w: 



(m le^^^l n> = exp (^-^) (.7)">~"< ^^T'"^^ (7') > 

n> = max (n, m) : 

^ ^ (6.131) 
n< = min (n, m) , 

7 = g-y/ 1^) quantum states of the one-dimensional oscillator with the frequency 

w, L^n^ [z) is the generalized Laguerre polynomial. 

System with mixed spins 

The correlation functions for a system with mixed spins can be explicitly derived by the 
generalization of Eqs. (I6.112p and (16.1291) to the case of the particles with the non-zero spin. 
We use the fact, that the derivation of Eqs. (I6.112p and (I6.129p . performed in this section. 
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does not depend on the concrete form of the energy spectrum Hence, in the formulae, 
derived above, the replacement should be made: 

|n) |n, a), iV^ Nn,a = <a«n,<T, (6.132) 

where a is the electron spin projection. Consequently, the matrix elements {(^^'^^)^^ are 
replaced by 

(^^^■*)n.n ^ ^ l^^^""! = (^^^■'')„.n " (6-133) 

Taking into account Eqs. ( I6.132p and 06.1331) . the two-point correlation function fl6.130p 
becomes 



n,m (Ti,(j2 

+ E E I (^^''■'') nm I ' " ^-)] (^n,. (l - Nm,.) ) • (6. 134) 



n,m a 



h 

The averages (^iVn o-) and ( A^n,o-i^m,(T2 ) are, respectively, one-particle and two-particle 
distribution functions, 



Nn,a)=fi{n,(7\N^,(3), (6.135) 
iVn,<xiVn',a') = /2 (n, a; n', a'| {iVj , (3) . (6.136) 



The one-electron distribution function fi (n, a\N„, /3) is the average number of electrons with 
the spin projection cr at the ra-th energy level, while the two-electron distribution function 
/2 (n, a; n', (j'l {N^} , (3) is the average product of the numbers of electrons with the spin 
projections a and a' at the levels n and n'. These functions are expressed through the 
following integrals [see ( 10271) . f l6:T28D ]: 

/i (n, a\N^, (5) = ^ / / (£„, 6) $ (9, (3, N^) dO, (6.137) 



f2{n,a-n',a'\{K},f3) 



— TT 

A (n, a|iV,, /3) /i (n, a'\N,.,P) , li a' ^ a 



(6.138) 

with the notations 



$ (^, /3, iV,) = exp 



5^1n(l + e^^+«-'^^")-iV.(e + 



.n=0 



(6.139) 
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1 



(6.140) 



exp {f3e - ^ - i6) + 1 

The function / {e, 6) formally coincides with the Fermi-Dirac distribution function of the 
energy e with the "chemical potential" (^ + i6) / f3. 

Here we consider the zero-temperature limit, for which the integrals (16.1371) and (I6.138P 
can be calculated analytically. The result for the one-electron distribution function is 



/i(n,a|/3,iV,)|^_ 



0, 

9h„ 



n < 

n > L„; 
n = L„. 



(6.141) 



According to (I6.14ip . is the number of the lowest open shell, and 

I i(n + l)(n + 2) {3D), 
[ n + 1 {2D) . 

is the degeneracy of the n-th shell. Nl^ is the number of electrons in all the closed shells 



with the spin projection a, 

.V,„ J^i'<^' + ^'<^'^^> (6.142) 
Si [ + {2D). 

The two-electron distribution function /2 {n, a; n', a' \ {Ncr} , /3) at T = takes the form 



f2{n,a;n',a'\P,{N.})\^^^ 

' /i (n, a\(3, N^) \^^^ h {n\ iV^O 1/3- 



0, 



n^^n' 01 a ^ a' 
(J = a' and n = n' < L(j\ 
u = o' and n = v! > L^j] 
a = a' and n = n' = L„. 



(6.143) 



In summary, we have obtained the following expression for ^ (q, t| {N^r} , (3): 
~g{ci,r\{N^},f3)= ^ (e''^-)^^ (e^'^-)^,^, (n, a; n', a'| {iVj , /3) 

n,(T,n',(T' 
n,n',(T 

X [A {n, a\ {iVj ,f3)-f2 {n, a; n', a\ {N^} , /?)] . 



— {Sn — Sn') 



(6.144) 



This formula is valid for both closed and open shells. The correlation functions derived in 
this subsection are used both for the calculation of the ground-state energy and, a shown 
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below, for the calculation of the optical conductivity of an A^-polaron system in a quantum 
dot. 



E. Many-polaron ground state in a quantum dot: extrapolatio n to the homoge- 
neous Umit and comparison to the results for a polaron gas in bulk jsJ] 



The correlation function given by Eq. 06.1441) can be subdivided as 



~g (q, r| {iVj , P) = ~g, (q, r| {iVj , /3) + ^2 (q, r| {iVj , /3) 



(6.145) 



with 



n,n',(T 

X [U (n, a\ {iVj , /3) - /2 (n, a; n', a| {iVj 



r 

Vh 



(6.146) 



^2 



In accordance with the subdivision fl6.145p of the correlation function, we subdivide the 
Coulomb and polaron contributions: 



Ec = + E^'^ 



C ' 



E. 



(6.148) 
(6.149) 



e-ph — ^i-ph + ^e-ph- 

We have numerically checked whether the polaron contribution per particle E\^2phl-^ 
tends to a finite value at — t- 00. In Figs. [14] and [HI we have plotted the polaron 
contributions E^f^^f^/N as a function of for a quantum dot in ZnO and in a polar medium 
with a = 5, ?7 = 0.3, respectively. ^ As seen from Fig. [TH the polaron contribution E^^}^y^/N 
in ZnO as a function of A^ oscillates taking expressed maxima for A^ corresponding to the 
closed shells A^ = 2, 8, 20, 40, . . .. There exist kinks of E^^}^f^/N at A^ corresponding to the 
half-filled shells, but these kinks are extremely small. In the case of the medium with a = 5, 
7] = 0.3, for r* = 20 (what corresponds to the density no ^ 1.14 x 10^® cm~^), the polaron 



Since, as discussed above, for a single polaron only the whole polaron contribution E(.-ph = E^^} + E^'^^ 



has a physical meaning, the plots for -Egj'p/j in Figs. 3 and 4 start from N = 2. The total polaron 
contribution Ee-ph ioi N = 1 is plotted below, in Fig. 9. 
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contribution e'^^}^^/N oscillates taking maximal values at the numbers of fermions, which 
correspond to the closed shells for a spin-polarized system with parallel spins. 

In Figs. [H] and [151 the dashed curves are the envelopes for local maxima (closed shells) 
and local minima of E^^}^^/ N IWe see that when these envelopes are extrapolated to larger 
number of fermions, the distance between the envelopes decreases. Therefore, the magnitude 
of the variations of e'^^}^^/N related to the shell filling diminishes with increasing A^, and 
it is safe to suppose that in the limit of large A^, the envelopes tend to one and the same 
value. That value corresponds to the homogeneous ("bulk") limit limjv-5.00 (-^e-p/i/^j- 
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= 0.849 
= 0.4908 



I Closed shells 
I Half-filled shells 
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FIG. 14: Polaron contribution E(,_^f^/N as a function of N for a ZnO quantum dot. The material 
parameters for ZnO are taken from Ref. [8{|]. The value r* = 2 corresponds to no = 4.34 x 10^^ 
cm~^. Inset: the total spin as a function of N . 
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FIG. 15: Polaron contribution E^^l^^/N 3jS Oil function of for a quantum dot of a polar medium 
with a = 5, 7/ = 0.3. The value r* = 20 corresponds to no = 1.14 x 10^^ cm^^. Inset: the total 
spin as a function of A^. 



In Fig. [ini we compare the polaron contribution E^^} ,/N calculated within our varia- 



tional path-integral method for different numbers of fermions with the polaron contrib ution 
to the ground-state energy per particle for a polaron gas in bulk, calculated (i) in R ef [101 1 
within an intermediate-coupling approach (the thin solid curve), (ii) in Ref. 1Q2|, using a 



151 



variational approach developed first in Ref. 89|]. As seen from this figure, our all-coupling 
variation al method provides lower values for the polaron contribution than those obtained 
in Refs. 



101 



102|. The d iffer ence between the polaron contribution calculated within our 



method and that of Ref. 



10l| is smaller at low densities and increases in magnitude with 



increasing density. The differ ence between the polaron contribution calculated withi n our 
method and that of Ref. 102l | very slightly depends on the density. The result of Ref. 102 1 



becomes closer to our result only at high densities. 
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FIG. 16: Polaron contribution to the polaron ground-state energy per particle E^^}^^^ in an A^- 
polaron quantum dot as a function of the effective density. The parameters are taken for ZnO (see 
Ref. [Q]): a = 0.849, Eq = 8.15, e^o = 4.0, /iwlo = 73.27 meV, rrn, = 0.24me, where nie is the 



electron mass in the vacuum. This polaron contribution is compared wit 
to the ground-state energy of a polaron gas in bulk calculated in Refs. 



1 the polaron contribution 



1 tri e yoV i 
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152 



F. Optical conductivity 



In Ref. 96| we have extended the memory-function approach to a system of arbitrary- 
couphng interacting polarons confined to a parabohc confinement potential The optical 
conductivity relates the current J (t) per electron to a time-dependent uniform electric field 
E (t) in the framework of linear response theory. Further on, the Fourier components of the 
electric field are denoted by E (w) : 

1 r°° 

E(t) = — / E{uj)e-"^'duj, (6.150) 

and the similar denotations are used for other time-dependent quantities. The electric 
current per electron J (t) is related to the mean electron coordinate response R (t) by 

dR (t) 



J (t) = -e- 



dt 



(6.151) 



and hence 

J (w) = iewR (w) . (6.152) 

Within the linear-response theory, both the electric current and the coordinate response are 
proportional to E (cu): 

a{uj). 



leixJ 



-E(w) 



(6.153) 



where a (w) is the conductivity per electron. Because we treat an isotropic electron-phonon 
system, a {u) is a scalar function. It is determined from the time evolution of the center-of- 
mass coordinate: 

= • (6.154) 



ii=i 



The symbol {{{*))) g denotes an average in the real-time representation for a system with 
action functional S: 



(((•)))< 



X X 

rfx j d^,j d^,j D^it) j Dx' (t) ei^t^W'^'W] (•) (xo|p(to)|x; 



(6.155) 

where (xq |p (to) I Xg) is the density matrix before the onset of the electric field in the infinite 
past (to — )■ — oo). The corresponding action functional is 

5[x(t),x'(t)] = 



Le [k (t) , X (t) , t) - Le (^X (t) , X' (t) , t j J dt' - im [X (t) , x' (t)] , 

(6.156) 
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where Le (x, x, t) is the Lagrangian of interacting electrons in a time- dependent uniform 
electric field E (t) 

/w.x2 2 . TV, 

(M. ^ E E - ^ - ■ E w - E E E ,.,|.:_.„| - 

a j=l ^ ^ a,a' j=l 1=1 °° ' ^''^ ''"^ ' 

The influence phase of the phonons 



(6.157) 



,2 ' 



['V - Pa [S 



$[x(.),x'(.)] = -5^^ J dsj ds'[p^ 

^ — oo — oo 

X [t^ - s') (s') - T^^ - is')] (6.158) 

describes both a retarded interaction between different electrons and a retarded self- 
interaction of each electron due to the elimination of the phonon coordinates. This functional 
contains the free-phonon Green's function: 

iujt —iuit 

(t) = — + — . (6.159) 

The equation of motion for R (t) is 

m^^^^l^ + ^b^oR W + eE (t) = Fp;, (t) , (6.160) 
where ¥ph it) is the average force due to the electron-phonon interaction, 

F,.(i) = -Re5^^^P j dsT:^Jt-s){{p^{t)p^^{s))),. (6.161) 

1 -oo 

The two-point correlation function ((pq (t) p_q (s)))5 should be calculated from Eq. fl6.155p 
using the exact action (16.1561) . but like for the free energy above, this path integral cannot 
be calculated analytically. Instead, we perform an approximate calculation, replacing the 
two-point correlation function in Eq. (I6.16ip by ((pq (t) p_q (s)))5o, where So [x (t) , x' (t)] 
is the action functional with the optimal values of the variational parameters for the model 
system considered in the previous section in the presence of the electric field E(t). The 
functional Sq [x (t) , x' (t)] is quadratic and describes a system of coupled harmonic oscillators 
in the uniform electric field E (t). This field enters the term — eE (t) ■ Ylf=i in the 
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Lagrangian, which only affects the center-of-mass coordinate. Hence, a shift of variables to 
the frame of reference with the origin at the center of mass 



X, (t) =5„(t) + R(t) 

x;(t) = i;(t) + R(t) 



(6.162) 



results in 



{{pAt) P-As))) So = ((p. WP-.(^)))5oU=o^^''^''^*^'''^^^^- (6.163) 
This result (16.1631) is valid for any quadratic model action 5*0. 

The applicability of the parabolic approximation is confirmed by the fact that a self- 
induced polaronic potential, created by the polarization cloud around an electron, is rather 
well described by a parabolic potential whose parameters are determined by a variational 
method. For weak coupling, our variational method is at least of the same accuracy as the 
perturbation theory, which results from our approach at a special choice of the variational 
parameters. For strong coupling, an interplay of the electron-phonon interaction and the 
Coulomb correlations within a confinement potential can lead to the assemblage of polarons 
in multi-polaron systems. Our choice of the model variational system is reasonable because 
of this trend, apparently occurring in a many-polaron system with arbitrary N for a finite 
confinement strength. 

The correlation function (pq (t) p_q (s))^^ corresponds to the model system in the ab- 
sence of an electric field. For t > s, this function is related to the imaginary-time correlation 
function Q (q, r| {N^} , (3) , described in the previous section: 

((pq it) p_q (5)))5Je=o,>. = ^ ^ - I {^4 , /3) ■ (6.164) 
Using the transformation ( 16.162^ and the relation ( 16.164^ . we obtain from Eq. ( 16.16ip 

2 * 

F,,(t) = -ReJ]^^P I T:^Jt-.)e^''-[^W-^W]^(q,^(t-.)|{iV4,/3)ci.. (6.165) 

^ -oo 

Within the linear- response theory, we expand the function e^'^ t^'^*^"^'^*^] in Eq. (I6.165P 
as a Taylor series in [R (t) — R (s)] up to the first-order term. The zeroth-order term gives 
no contribution into Fph (t) due to the symmetry of \Vqf and of /q {t — s) with respect to 
the inversion q — q. In this approach, the Cartesian coordinates of the force (j = 1, 2, 3) 
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become 

t 



X Im [T^^ {t-s)g (q, t{t-s)\ {iVj , /?)] rfs. (6.166) 
Further on, we perform the Fourier expansion: 

1 r°° 

R(t) = — / RiLo)e-"^*'du. (6.167) 
In Eq. fl6.166p . we make the replacement 

T = t — S, =^ S = t — T, 

what gives 

3 2 °° 

(Fp/. (^)). = E E ^^^^ / dr [R, (t) -R,{t- r)] Im [T^^^ (r) ^ (q, zt\ {iVj , /?)] 

k=l q 

3 2 °° 

= E E /°° (-) / * (1 - K„ (-) /. 

fc=i q "^-"^ ;^ 

1 /""^ 

= - j JcoF,iuj)e-^-\ 
where the Fourier component of the force is 

3 2 °° 

(Fp. (^)), = E E ^^^^^ I rft (1 - e^-) Im [T^^ (r) ^ (q, ^t\ {iVj , /?)] R, 



/ql^ 



A:=l q 



A: [UJ) ■ 

(6.168) 



The expression f l6.168p can be written down as 

3 

(Fph iuj))j = -rrib Xjk (uj) Rk (uj) , (6.169) 

k=l 

where Xjk {^) sue components of the tensor 

2 °° 

(^) = E ^^i^ / (^'"* - 1) [^^.o (t) G (q, ^^1 {N.} , /?)] • (6.170) 

In the abstract tensor form, Eq. fl6.169p is 

¥{uj) = -'^{uj)K{uj). (6.171) 
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In particular, for the isotropic electron-phonon interaction and in the absence of the magnetic 
field, the tensor ^ (w) is proportional to the unity tensor I, 



^(a;)=xMl, (6.172) 



where x {^) is the scalar memory function: 

|2 „2 



1 



Let us perform the Fourier transformation of the equation of motion f l6.160p : 

nib i^l - oj^) R (oj) + eE (u) = Fpt (u) . (6.174) 
With Eq. f l6.17ip . this equation takes the form 



rrib {^l - w^) R (io) + eE (w) = -m^ Y (oo) R 



[UJ 



4 

nib [co^ -nl-^ (tu)] R (to) = eE (w) . (6.175) 
Comparing Eqs. fl6.153p and f l6.175p between each other, we find that 



nib \uj^ - - V i^)] (io) = eE (u) , 

ieu) 



so that 



nib 



In the case when Eq. (16.1720 is valid, we obtain the conductivity in the scalar form 

nib [w^ -il^-x Ml 

The real part of the conductivity is 

ie'^u [(w^ — fig) — Re X (c^) + ^ Ini X (w)] 



Re cr (w) = Re 

I [{u^ - Ql) -Rex (co)? + [Im x {uj)f 
e^uj Imx(w) 
= ~^b [(a;2 _ ^2) _ Rex (u;)]' + [Imx (a;)]'' 
In summary, the optical conductivity can be expressed in terms of the memory function 



X{uj) (cf. Ref. |48|), 

e^ w Im X {oj 

m,[uj^-ai-Rex{uj)f + [Imx (a;)] 
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where x i!^) is given by 

1 ^ 

It is worth noting that the optical conductivity fl6.176p differs from that for a translationally 
invariant polaron system both by the exphcit form of x {'^) ^"^^ by the presence of the term 
in the denominator. For a — )■ 0, the optical conductivity tends to a 5-like peak at uj = VLq, 

limRea(u;) = 5 iw - VLq) . (6.178) 

For a translationally invariant system i7o — 0, and this weak- coupling expression fl6.178p 
reproduces the "central peak" of the polaron optical conductivity 



The further simplification of the memory function (16.1771) is performed in the following 
way. With the Frohlich amplitudes of the electron-phonon interaction, we transform the 
summation over q to the integral and use the Feynman units {h= 1, wlo = 1? '^fe = 1); in 
which |Vq|^ = We also use the fact that in an isotropic crystal, Q (q, it\ {N^j} , /3) = 

Q (g, it\ {N„} , (3). As a result, we find 

oo 

X I (e^-* - 1) Im [T^^ (t) Q (g, it\ {iVj , /?)] 



oo 

= 1^ \^ q'dq I dt (e^-* - 1) Im [T^^^ (t) Q (g, zt| {iV^ , (3)] . 



In the zero-temperature case, T* (t) — e~**, and we arrive at the expression 



oo 

g^rfg J dt (e*"* - 1) Im [e~''g (g, it| {A^,} , (3)] . (6.179) 



Substituting the two-point correlation function Q {q,it\ {Ncr} , P) with the one-electron 
(16.1411) and the two-electron (I6.143P distribution functions into the memory function (I6.179P 
and expanding Q (g, it\ {N^} , /3) in powers of e~™*, e"**^^* and e~*^^*, the integrations over g 
and t in Eq. (I6.179P are performed analytically. The similar transformations are performed 
also in the 2D case. As a result, the memory function (16.177^ is represented in the unified 
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form for 3D and 2D interacting polarons: 

^ ^ ' e^+O SttA^ 14/ \ A ) 



X 



oo oo oo 



-1) 



P3 



P1=0 P2=0 P3=0 
OO oo 



at 



pi 



P2 



P3 



X 



E E E ^1 {^-> ' - /2 (r^, a; m, a| {iVj , /?)] |^^^ 



m=0 n=0 <T 



m n 



-1) 



^^LO +Pi f^i +P2 f^2 + (P3 -rri+n)to 
n— m+i+fc 



Z=0 fc=?i— m+i 

n — k J \k — I — n + m 



T {pi+P2+P3 + k + l + l) 

kill [wAJ 



l+k 



+ 



1 



1 



CJ-WLO--{Plf^l+P2f22+P3«')+i£ UJ+LJi^o+piQ,\+p20,2+P3'W+i£ 

+v 



X 



<^LO+Plf2l+P2f22+P3W 

E E E (n, a; m, a'\ {N,} , 



m=0 n=0 CTjct' 
71 m / -r \k+l 



fc=o ;=o 
n + D - l\fm + D - I 
n — k I \ m — I 



(6.180) 



where D = 2, 3 is the dimensionahty of the space, V denotes the principal value, A is defined 
as A = [XlLi c^l/^i + {N — 1) /w] /N, Qi,Q2, and w are the eigenfrequencies of the model 
system, ai and 02 are the coefficients of the canonical transformation which diagonalizes the 
model Lagrangian fl6.32p . 



1. Selected results: the manifestations of the shell filling in optical conductivity 

The shell filling schemes for an A^-polaron system in a quantum dot can manifest them- 
selves in the spectra of the optical conductivity. In Fig.[T71 optical conductivity spectra for 
= 20 polarons are presented for a quantum dot with the parameters of CdSe: a = 0.46, 
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Q 



Tj = 0.656 [5[ and with different values of the confinement energy hflQ. ^ In this case, 
the spin-polarized ground state changes to the ground state satisfying Hund's rule with 
increasing hUo in the interval 0M21H* < < 0M22H*. 



o 
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FIG. 17: Optical conductivity spectra of A'^ = 20 interacting polarons in CdSe quantum dots with 
a = 0.46, r] = 0.656 for different confinement energies close to the transition from a spin-polarized 
ground state to a ground state obeying Hund s rule. Inset: the first frequency moment (cj) of the 
optical conductivity as a function of the confinement energy. (From Ref. [96].) 



In the inset to Fig.HTJ the first frequency moment of the optical conductivity 

jp Kea[uj)duj 

as a function of MIq shows a discontinuity, at the value of the confinement energy corre- 
sponding to the change of the shell filling schemes from the spin-polarized ground state to 
the ground state obeying Hund's rule. This discontinuity might be observable in optical 
measurements. 



* For the numerical calculations, we use effective atomic units, where /i, the electron band mass and 
e/ j£ao have the numerical value of 1. This means that the unit of length is the effective Bohr radius 



i|j = f^Eaal (^Tifce^), while the unit of energy is the effective Hartree iJ* — m^^e^ j {T?e 



oo J 
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The shell structure for a system of interacting polarons in a quantum dot is clearly 
revealed when analyzing the addition energy and the first frequency moment of the optical 
conductivity in parallel. In Fig [TSl we show both the function 



(6.182) 



and the addition energy 



A (iV) = E°{N + 1)- 2E^ (N) + E°{N - 1) 



(6.183) 



for interacting polarons in a 3D CdSe quantum dot. 
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FIG. 18: The function B (N) and the addition energy A (N) for systems of interacting polarons in 



CdSe quantum dots with a = 0.46, r] = 0.656 for = O.IH* . (From Ref. 



As seen from Fig [18], distinct peaks appear in (A^) and A (A^) at the "magic numbers" 
corresponding to closed-shell configurations at = 8, 20 and to half-filled-shell configura- 
tions at = 5, 14. We see that each of the peaks of G (A^) corresponds to a peak of the 
addition energy. The filling patterns for a many-polaron system in a quantum dot can be 
therefore determined from the analysis of the first moment of the optical absorption for 
different numbers of polarons. 
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VII. VARIATIONAL PATH-INTEGRAL TREATMENT OF A TRANSLATION 
INVARIANT iV-POLARON SYSTEM 



A. The many-polaron system 

In the present section, the ground-state properties of a translation invariant A^-polaron 



system are theoretically studied in the framework o: 
'or identical particles, using a further development 



96, 



98 



l3- 



the variational path-integral method 



103( 1 of the model introduced in Refs. 



In order to describe a many-polaron system, we start from the translation invariant A^- 
polaron Hamiltonian 

N 2 N N 2 / ^ 



^ = E If: + ^ E E r^TH + E ^-oaia, + E E v^--^''' + h-c , (r.i) 

ZjIIL z, too I'-i , \ • 1 1 



where m is the band mass, e is the electron charge, culo is the longitudinal optical (LO) 
phonon frequency, and Vk are the amplitudes of the Frohlich electron-LO-phonon interaction 

^ k \ V ) \2muji^o) ' 2nuJi,o \ ^ } \^oo eo/ ' 

with the electron-phonon coupling constant a > 0, the high-frequency dielectric constant 
eoo > and the static dielectric constant eo > 0, and consequently 

!t^^n('^\'\^aV-2<{^\'\u. (7.3) 

eoo \ m J \nwLoJ 

In the expression (17. 3p . H* is the effective Hartree 

H* = ^, a*B = ^^ (7.4) 

where a*^ is the effective Bohr radius. The partition function of the system can be expressed 
as a path integral over all electron and phonon coordinates. The path integral over the 



phonon variables can be calculated analytically [4l|. Feynman's phonon elimination tech- 
nique for this system is well known and leads to the partition function, which is a path 
integral over the electron coordinates only: 
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where r = {ri, ■ ■ ■ , r^} denotes the set of electron coordinates, and §'Dr denotes the path 
integral over all the electron coordinates, integrated over equal initial and final points, i.e. 

.r(/3)=r 



'f(0)=r 

Throughout this paper, imaginary time variables are used. The effective action for the 
A^-polaron system is retarded and given by 

^ I m / dYj (r) \ ^ 1 
\2^V dr ) ^2^L.e^\v,{r) -n (r) 



j=i ^ ^ j=i i=i,^j 

+ lf fflY. |H^-'-(-.W-.'"> '°^"'^." - I" - i.dr. (7.6) 

Note that the electrons are fermions. Therefore the path integral for the electrons with par- 
allel spin has to be interpreted as the required antisymmetric projection of the propagators 
for distinguishable particles. 

We below use units in which ^ = 1, m = 1, and wlo = 1- The units of distance and 

1 /2 

energy are thus the effective polaron radius [h/ (mc^Lo)] aiid the LO-phonon energy hui\,o. 
B. Variational principle 



For distinguishable particles, it is well known that the Jensen- Feynman inequality [41 1 
provides a lower bound on the partition function Z (and consequently an upper bound on 
the free energy F) 

'(7.7) 
I — 

for a system with real action 5* and a real trial action 5*0. The many-body extension (Ref. 971 . 



105!]) of the Jensen- Feynman inequality, requires that the potentials are symmetric with 
respect to all particle permutations, and that the exact propagator as well as the model 
propagator are defined on the same state space. Within this interpretation we consider the 
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following generalization of Feynman's trial action 



with the variational frequency parameters v,w,u. 

Using the explicit forms of the exact (I7.6P and the trial (I7.9p actions, the variational 
inequality (17. 8p takes the form 



U / ^ 1 



dr 

ry.: 





2 I 2 1 i>B I ^ 



N(3 

2/3 70 io V /, sinh>,o 

(7.10) 

In the zero-temperature limit (/3 — )■ oo), we arrive at the following upper bound for the 
ground-state energy {N^, N^) of a translation invariant N-polaTon system 

with 

Q / \ 2 Q 1 -I 

E^ar (iVt, N^\v, w, u) = ^ ^""'J^' - + W + (iVj 

+ Ecu W + Ec\\ {N^) + Ecu (A^t, N^) 

+ (iVt) + (iVj + (iV^, iVj , (7.11) 

where Ej? (A^) is the energy of spin-polarized fermions confined to a parabolic poten- 
tial with the confinement frequency Ec\\ [^N^(^i)) is the Coulomb energy of the electrons 
with parallel spins, Ec^^ {N^, N^) is the Coulomb energy of the electrons with opposite 
spins, Ea\\ (^t(4)) is electron-phonon energy of the electrons with parallel spins, and 
Ea-\i (iV|, N\^) is the electron-phonon energy of the electrons with opposite spins. 
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C. Results 



Here, we discuss some results of the numerical minimization of E^jar {N^, N^lv, w, u) with 
respect to the three variational parameters v, w, and u. The Frohlich constant a and the 
Coulomb parameter 

cto = —1= = with - = — (7.12) 

v2 



— with — = — ^ 



characterize the strength of the electron-phonon and of the Coulomb interaction, obeying 
the physical condition a > [see (17.31) ]. The optimal values of the variational parameters 
f,w,and u are denoted Vop,Wop,Biid Uop, respectively. The optimal value of the total spin 
was always determined by choosing the combination {N^, N^) for fixed N = + iVj, which 
corresponds to the lowest value (N) of the variational functional 



opi '^opi 



1 l^op) . 



(7.13) 




FIG. 19: The "phase diagrams" of a translation invariant A^-polaron system. The grey area is 
the non-physical region, for which a > oq. The stab ility region for each number of electrons is 



determined by the equation ac < ca < a^. (From Ref. 



103]. 



In Fig. Uni the "phase diagrams" analogous to the bipolaron "phase diagram" of Ref. 



106| are plotted for an A^-polaron system in bulk with = 2, 3, 5, and 10. The area where 
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ao < a is the non-physical region. For a > ao, each sector between a curve corresponding 
to a well defined N and the line indicating ao = a shows the stability region where Uop 7^ 0, 
while the white area corresponds to the regime with Uop = 0. When comparing th e sta bility 
region for = 2 from Fig. with the bipolaron "phase diagram" of Ref. 106| . the 
stabil ity r egion in the present work starts from the value ac ~ 4.1 (instead of etc ~ 6.9 in 



Ref. 



106| ). The width of the stability region within the p resen t model is also larger than 



the width of the stability region within the model of Ref. 



106|. Also, the absolute values 



of the ground-state energy of a two-polar on sy stem given by the present model are smaller 



than those given by the approach of Ref. 



106 |. 



The difference between the numerical results of the present work and of Ref. lOq is 
due to the following distinction between the used model systems. The model system of Ref. 
106| consists of two electrons interacting with two fictitious particles and with eac h ot her 



through quadratic interactions. But the trial Hamiltonian given by Eq. (6) of Ref. 106| is 



not symmetric with respect to the permutation of the electrons. It is only symmetric under 
the permutation of the pairs "electron + fictitious particle". As a consequence, this trial 
system is only applicable if t he e lectrons are distinguishable, i.e. have opposite spin. In 



contrast to the model of Ref. 



106|, the model used in the present paper is described by the 



trial action (9), which is fully symmetric with respect to the permutations of the electrons, 
as is required to describe identical particles. 

The "p hase diagrams" for N > 2 deinonstrate the existence of stable multipolaron states 
(see Ref. jl07 |). As distinct from Ref. 107 |. here the ground state of an A^-polaron system 
is investigated supposing that the electrons are fermions. As seen from Fig. [201 for N > 2, 
the stability region for a multipolaron state is narrower than the stability region for N = 2, 
and its width decreases with increasing A^. 

A consequence of the Fermi statistics is the dependence of the polaron characteristics 
and of the total spin of an A^-polaron system on the parameters {a, ao,N). In Fig. [20l we 
present the ground-state energy per particle, the confinement frequency Uop and the total 
spin S" as a function of the coupling constant a for ao/a = 1.05 and for a different numbers 
of polarons. The ground-state energy turns out to be a continuous function of a, while 
Wopand S reveal jumps. For N = 2 (the case of a bipolaron), we see from Fig. [20] that the 
ground state has a total spin S = for all values of a, i. e., the ground state of a bipolaron 
is a singlet. This result is in agreement with earlier investigations on the large-bipolaron 
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FIG. 20: The ground-state energy per particle (a), the opthnal value Uop of the confinement 
frequency (b), and the total spin (c) of a translation invariant A^-polaron system as a function of 
the coupling strength a for ao/a = 0.5. The vertical dashed lines in the panel c indicate the critical 
values Oc separating the regimes of a > Oc, w here the multipolaron ground state with ujop ^ 



exists, and a < Oc, where uJop = 0. (From Ref. 



1031 1 . 



problem (see, e. g., [l08|). 

In summary, using the extension of the Jensen-Feynman variational principle to the 
systems of identical particles, we have derived a rigorous upper bound for the free energy of 
a translation invariant system of interacting polarons. The developed approach is valid 
for an arbitrary coupling strength. The resulting ground-state energy is obtained taking 
into account the Fermi statistics of electrons. 
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VIII. RIPPLONIC POLARONS IN MULTIELECTRON BUBBLES 



A. Ripplon-phonon modes of a MEB 

Spherical shells of charged particles appear in a variety of physical systems, such as 
fuUerenes, metallic nanoshells, charged droplets and neutron stars. A particularly interesting 
physical realization of the spherical electron gas is found in multielectron bubbles (MEBs) 
in liquid helium-4. These MEBs are 0.1 fim - 100 fim sized cavities inside liquid helium, 
that contain helium vapor a t vap or pressure and a nanometer-thick electron layer anchored 



to the surface of the bubble 



110l | . They exist as a result of equilibr ium between the surface 



tension of liquid helium and the Coulomb repu 



sion of the electrons 



lll| . Recently proposed 



experimental schemes to stabilize MEBs 112| have stimulated theoretical investigation of 
their properties. 

We describe the dynamical modes of an MEB by considering the motion of the helium 
surface ("ripplons") and the vibrational modes of the electrons toget her. In particular, we 
analyze the case when the ripplopolarons form a Wigner lattice 109 1 . 

First, we derive the Lagrangian of interacting ripplons and phonons within a continuum 
approach. The shape of the surface of a bubble is described by the function R {6, ip) = 

+ u {6,(f) , where u {6,!f) is the deformation of the surface from a sphere with radius 
i?b- The deformation can be expanded in a series of spherical harmonics Yim {0, ip) with 
amplitudes Qim, 

oo I 
1=1 m=—l 

We suppose that the amplitudes are small in such a way that ^/UT^^T) \Qim\ ^ -Rb- 

The ripplon contribution (Tf) to the kinetic energy of an MEB, and the contributions 
to the potential e nergy due to the surface tension (f/o-) and due to the pressure {Uy) were 
described in Ref. |ll3| : 

I 



1=1 m=—l 



= AnaRl + ^ (/^ + / + 2) |Q,„| 



An 



1=1 m=—l 
oo I 



(=1 m=—l 
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Here p ~ 145 kg/m^ is the density of liquid helium, a ~ 3.6 x 10"'* J/m^ is its surface 
tension, and p is the difference of pressures outside and inside the bubble. 

Expanding the surface electron density n {6, (f) in a series of spherical harmonics with 
amplitudes nim, 

oo I 

n (^, V;) = 5^ 5^ n^^y^^ (^, (^) , (8.3) 

/=0 m=-l 

the kinetic energy of the motion of electrons can be written as 

oo I 



T 



47rme-Rb . 2 

I ^Im I ; 



4) 



where me is the bare electron mass and is the number of electrons. Finally, the electrostatic 
energy [Uq) of the deformed MEB with a non-uniform surface electron density fl8.3p is 
calculated using the Maxwell equations and the electrostatic boundary conditions at the 
surface. The result is: 



p2 /V2 00 I 

1=1 m=—l 
2Ar2 °° ' /2 
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b l=\ m=-l ' V ' V 



EE 
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00 / 



/ + i 



-nimQ 
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1=1 m=—l 

with the dielectric constant of liquid helium e ^ 1.0572. The last term in Eq. (18. 5p describes 
the ripplon-phonon mixing. Only ripplon and phonon modes which have the same angular 
momentum couple to each other. After the diagonalization of the Lagrangian of this ripplon- 
phonon system, we arrive at the eigenfrequencies: 



ul (/) + (/) 



1/2 



(8.6) 



where u^- (/) is the bare ripplon frequency. 



Ur (/) 



' A^eRl / + £(/ + 1) 



1/2 



(8.7) 
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while Up (/) is the bare phonon frequency, 



(8.8) 



meRl / + £:(/ + !_ 
and 7 (/) describes the ripplon-phonon coupling: 

fig [imn^pKiJ l + e{l + iy ^ ' 

B. Electron-ripplon interaction in the MEB 

The interaction energy between the ripplons and the electrons in the multielectron bubble 
can be derived from the following considerations: (i) the distance between the layer electrons 
and the helium surface is fixed (the electrons find themselves confined to an effectively 2D 
surface anchored to the helium surface) and (ii) the electrons are subjected to a force field, 
arising from the electric field of the other electrons. For a spherical bubble, this electric field 
lies along the radial direction and equals 



E = -^er. (8.10) 

A bubble shape oscillation will displace the layer of electrons anchored to the surface. The 
interaction energy which arises from this, equals the displacement of the electrons times the 
force eE acting on them. Thus, we get for the interaction Hamiltonian 

= ^e|E| XM(fi,). (8.11) 
j 

Here u{Q) is the radial displacement of the surface in the direction given by the spherical 
angle fl; and Clj is the (angular) position operator for electron j. The displacement can be 
rewritten using (18. ip and we find 

H^nt = J]e|E| J]g,„rf™(fi,). (8.12) 

j £,m, 

Using the relation 

Qi,n. = + K-J^ (8-13) 

the interaction Hamiltonian can be written in the suggestive form 

H^nt = J2Y1 Me,mYe,U^j)ide^m + a^.J, (8.14) 

i,m j 
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with the electron-ripplon couphng amphtude for a MEB given by 



M, 



i,in 



(_l)("^-|™l)/2 



Ne^ h{i + l) 
2if V 2pRlu, 



U5) 



C. Locally flat approximation 



Substituting M^^^ into (18.141) . we get 

Hint 



/n(^+_i) 



(8.16) 



X 



_l\{m-\m\)/2 



Rh 



,m + Q-i jyj- 



In this expression, we consider the hmit of a bubble so large that the surface becomes flat 
on all length scales of interest. Hence we let R^ oo but keep i/Rb = q a constant. This 
means we have to let ^ — )■ oo as well. In this limit, 



hm Ye,o{9) = ^== 



sin[(£+l/2)^ + 7r/4], 



(8.17) 



and Yifl{9) varies locally as a plane wave with wave vector q = i/Rb- The wave function 
Yi^mi^j) / Rb is furthermore normalized with respect to integration over the surface (with 
total area AnRl). Thus, we get in the locally flat approximation 



int 



^ ^ Ne^ I hq 



(8.18) 



or 



q J 



hq 



Ma = e E 



(8.19) 



This corresponds in the limit of large bubbles to the interaction Hamiltonian expected for a 
fiat surface. 



D. Ripplopolaron in a Wigner lattice: the mean-field approach 



In their treatment of the electron Wigner lattice embedded i n a p olarizable medium such 
as a semiconductors or an ionic solid, Fratini and Quemerais 114| described the effect of 
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the electrons on a particular electron through a mean-field lattice potential. The (classical) 
lattice potential Viat is obtained by approximating all the electrons acting on one particular 
electron by a homogenous charge density in which a hole is punched out; this hole is centered 
in the lattice point of the particular electron under investigation and has a radius given by 
the lattice distance d. 

Within this particular mean-field approximation, the lattice potential can be calculated 
from classical electrostatics and we find that for a 2D electron gas it can be expressed in 
terms of the elliptic functions of first and second kind, E (x) and K (x), 

2e2 



Viat r) 



!^\d-r\E 


Ard 






{d-rf_ 



nd'^ 

+ {d + r) sgn [d — r) K 



Ard 



(d-rf 



1.20) 



Here, r is the position vector measured from the lattice position. We can expand this 
potential around the origin to find the small-amplitude oscillation frequency of the electron 
lattice: 

lim V^at (r) = + lrn,u;l,r' + O (r^) , (8.21) 

with the confinement frequency 



^iat=\l^^. (8.22) 

In the mean-field approximation, the Hamiltonian for a ripplopolaron in a lattice on a locally 
flat helium surface is given by 

H = + Viat (f ) + ^ f\hj{q)a^d^ 

q 

+ ^M,e-*'i-^(aq + a+J, (8.23) 
q 

where f is the electron position operator. 

Now that the lattice potential has been introduced, we can move on and include effects 
of the bubble geometry. If we restrict our treatment to the case of large bubbles (with 

> 10^ electrons), then both the ripplopolaron radius and the inter-electron distance d are 
much smaller than the radius of the bubble R^. This gives us ground to use the locally fiat 
approximation using the auxiliary model of a ripplonic polaron in a planar system described 
by ( ]8.23p . but with a modified ripplon dispersion relation and an modified pressing field. 
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We find for tlie modified ripplon dispersion relation in the MEB: 



^(q) = X-q' + ^q, (8.24) 

V p pRb 

where Rh is the equihbrium bubble radius which depends on the pressure and the number of 
electrons. The bubble radius is found by balancing the surface tension and the pressure with 
the Coulomb repulsion. The modified electron-ripplon interaction amplitude in an MEB is 
given by 

The effective electric pressing field pushing the electrons against the helium surface and 
determining the strength of the electron-ripplon interaction is 

E = -^er. (8.26) 
E. Ripplopolaron Wigner lattice at finite temperature 

To study the ripplopolaron Wigner lattice at finite temperature and for any value of the 



4l| . This varia- 



electron-ripplon coupling, we use the variational path-integral approach 
tional principle distinguishes itself from Rayleigh-Ritz variation in that it uses a trial action 
functional Striai instead of a trial wave function. 

The action functional of the system described by Hamiltonian (18. 23 p . becomes, after 
elimination of the ripplon degrees of freedom, 



2 







dr / daG^Ur - a)e'^-^'^^^-'^''^\ (8.27) 





with 

G„(r-.) = S^!!Mtz^^_W2W. (8,28) 

In preparation of its customary use in the Jensen-Feynman inequality, the action functional 
(I8.27P is written in imaginary time t = ir with /3 = l/^ksT) where Tis the temperature. 
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We introduce a quadratic trial action of the form 

h(5 



Strial — I dr 





-^jdrj daGUr - a)r(r) • r(a). (8.29) 



where M, w, and Q are the variationally adjustable parameters. This trial action corresponds 
to the Lagrangian 

^ .2 _ ^^2 _ K _ 

u 2 2 2 2^^ ^ ^ 

from which the degrees of freedom associated with R have been integrated out. This La- 
grangian can be interpreted as describing an electron with mass me at position r, coupled 
through a spring with spring constant k, to its lattice site, and to which a fictitious mass M 
at position R has been attached with another spring, with spring constant K. The relation 
between the spring constants in fl8.30p and the variational parameters if, f2 is given by 



w 



VK/m„ (8.31) 



n = ^/{K + K)/m^. (8.32) 

Based on the trial action Stnaii Feynman's variational method allows one to obtain an 
upper bound for the free energy F of the system (at temperature T) described by the action 
functional 5* by minimizing the following function: 

F = Fq — — {S — Strial) 5 (8.33) 

with respect to the variational parameters of the trial action. In this expression, Fq is the 
free energy of the trial system characterized by the Lagrangian Cq, P = l/{ki,T) is the 
inverse temperature, and the expectation value {S — Striai) is to be taken with respect to 
the ground state of this trial system. The evaluation of expression f l8.33p is straightforward 
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though lengthy. We find 



/3 



2 sinh 
2 sinh 
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2 sinh 



^ ^ a^^^^i coth ( 
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X exp 
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sinh(Snj/3/2) 

In this expression, Jq and Ji are Bessel functions of imaginary argument, and 



1 1 bp i £0 



^2 



Finally, Qi and ^2 are the eigenfrequencies of the trial system, given by 
^1,2 = o 



(8.34) 



(8.35) 
(8.36) 

(8.37) 



Optimal values of the variational parameters are determined by the numerical minimization 
of the variational functional F as given by expression fl8.34p . 



F. Melting of the ripplopolaron Wigner lattice 



The Lindemann melting criterion 115| states in general that a crystal lattice of objects 
(be it atoms, molecules, electrons, or ripplopolarons) will melt when the average motion 
of the objects around their lattice site is larger than a critical fraction 6q of the lattice 
parameter d. It would be a strenuous task to calculate from first principles the exact value 
of the critical fraction 60, but for the particular case of electrons on a helium surface, we can 
make use of an experimental determination. Grimes and Adams 116| found that the Wigner 
lattice melts when F = 137±15, where F is the ratio of potential energy to the kinetic energy 
per electron. At temperature T the average kinetic energy in a lattice potential Viat is 



kin 



— - — coth 



2kBT 



(8.38) 
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and the average distance that an electron moves out of the lattice site is determined by 

From this we find that for the melting transition in Grimes and Adams' experiment 116| . 
the critical fraction equals 5q_~0.13. This estimate is in agreement with previous (empirical) 



estimates yielding 60 ~ 0.1 117l |. and we shall use it in the rest of this section. 



Within the approach of Fratini and Quemerais 114l |. the Wigner lattice of (rip- 
plo)polarons melts when at least one of the two following Lindemann criteria are met: 



6r = ^^^'"^^ > 5o, (8.40) 



= > ^0- (8.41) 



d 

where p and Rcms are, respectively, the relative coordinate and the center of mass coordinate 
of the model system fl8.30p : if r is the electron coordinate and R is the position coordinate 
of the fictitious ripplon mass M, this is 

^ mer + MR ^ 

^cms = , ; p = r - R. 8.42 

The appearance of two Lindemann criteria takes into account the composite nature of (rip- 
plo)polarons. As follows from the physical sense of the coordinates p and Rcms, the first 
criterion f l8.40p is related to the melting of the ripplopolaron Wigner lattice towards a rip- 
plopolaron liquid, where the ripplopolarons move as a whole, the electron together with 
its dimple. The second criterion f l8.4ip is related to the dissociation of ripplopolarons: the 
electrons shed their dimple. 

The path-integral variational formalism allows us to calculate the expectation values 
(Rims) and (p^) with respect to the ground state of the variationally optimal model system. 
We find 



/r)2 \ 

\^CTns / 



hw 



4 



^Q\{v? - VlI) coth(/i(]2/3/2)/fi2] , (8.43) 
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X [nl{w^ - nl) coth 

- w'^) coth(nn2^/2)] . 



(8.44) 



Numerical calculation shows that for ripplopolarons in an MEB the inequality i^i ^ w is 
fulfilled {w/Qi ~ 10^'^ to 10^^) so that the strong-coupling regime is realized. Owing to this 
inequality, we find from Eqs. fl8.43p . fl8.44p that 

(RL.) « (p'> . (8.45) 

So, the destruction of the ripplopolaron Wigner lattice in an MEB occurs through the 
dissociation of ripplopolarons, since the second criterion fl8.4ip will be fulfilled before the 
first fl8.40p . The results for the melting of the ripplopolaron Wigner lattice are summarized 
in the phase diagram shown in Fig. [211 
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FIG. 21: The phase diagram for the spherical 2D layer of electrons in the MEB. Above a critical 
pressure, a ripplopolaron solid (a Wigner lattice of electrons with dimples in the helium surface 
underneath them) is formed. Below the critical pressure, the rip plop olaron solid melts into an 
electron liquid through dissociation of ripplopolarons. (From Ref. 1091 ].) 
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For every value of A^, pressure p and temperature T in an experimentally accessible 
range, this figure shows whether the ripplopolaron Wigner lattice is present (points above 
the surface) or molten (points below the surface). Below a critical pressure (on the order 
of 10^ Pa) the ripplopolaron solid will melt into an electron liquid. This critical pressure 
is nearly independent of the number of electrons (except for the smallest bubbles) and is 
weakly temperature dependent, up to the helium critical temperature 5.2 K. This can be 
understood since the typical lattice potential well in which the ripplopolaron resides has 
frequencies of the order of THz or larger, which correspond to ~ 10 K. 

The new phase that we predict, the ripplopolaron Wigner lattice, will not be present for 
electrons on a flat helium surface. At the values of the pressing field necessary to obtain a 
strong enough electron-ripplon coupling, the flat helium surface is no longer stable against 
long- wavelength deformations jlisl. Multielectron bubbles, with their different ripplon dis- 
jers ion and the presence of stabilizing factors such as the energy barrier against fissioning 
119| . allow for much larger electric fields pressing the electrons against the helium surface. 
The regime of A^, p, T parameters suitable for the creation of a ripplopolaron Wigner lattice 
lies within the regime that would be a chievable in recently proposed experiments aimed at 
stabilizing multielectron bubbles |ll2]. The ripplopolar on Wigne r lattice and its melting 



transition might be detected by spectro scopi c techniques |ll6l . ll2Cll | probing for example the 



transverse phonon modes of the lattice 
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